Constant size and mobility of grains in the Oligo-Miocene megafan deposits in the Swiss Molasse Basin despite long-term changes in tectonic processes, climatic conditions, and sediment flux Open Access

In a foreland basin, the stratigraphic archives of alluvial fans made up of coarse-grained (> 2 mm) fluvial material have been used to infer the link between the processes operating in the orogen and the sedimentation dynamics in the basin and on the fans themselves. Among the various parameters, the size of grains in conglomerates and their grain-size distributions (GSDs) have received much attention where the aim has been to extract information on the sediment transport dynamics and controls thereof from the stratigraphic record (Paola et al. 1992; Duller et al. 2010; Armitage et al. 2011; Whittaker et al. 2011; Schlunegger and Castelltort 2016; Sharma et al. 2024). Core concepts for such interpretations of alluvial-fan stratigraphy and for reconstructing sediment transport dynamics are based on grain size (Paola and Mohrig 1996; Duller et al. 2010; Litty et al. 2016; Schlunegger and Garefalakis 2018), and they arise from the results of flume experiments (e.g., Meyer-Peter and Müller 1948; Wong and Parker 2006) and numerical, as well as conceptual models (e.g., Parker 1991; Paola et al. 1992; Whipple et al. 1998; Fedele and Paola 2007; Armitage et al. 2011). Based on these works, research has converged to a view that allogenic, or external, mechanisms such as changes in tectonic processes, shifts in paleoclimate, and modifications of exposed bedrock with different erodibility influence the sediment dynamics, including: the efficiency to generate sediment in the orogen, the supply of material to the foreland basin, and the formation of accommodation space in the basin itself (e.g., Flemings and Jordan 1989; Jerolmack and Paola 2010; Armitage et al. 2011; Whittaker 2012; Hajek and Straub 2017). Yet, the sensitivity of alluvial-fan stratigraphies and particularly of grain size to these driving forces has remained a topic of debate (Dorn 2009; Dade and Verdeyen 2007; Ventra and Nichols 2014; Meek et al. 2020; Griffin et al. 2023). For instance, models show that higher uplift rates in the source area, which is an allogenic driving force, result in an increase in surface erosion rates and a larger sediment flux to the foreland basin (Flemings and Jordan 1989; Paola et al. 1992). Yet, such a flux signal should also be associated with the supply of material with larger grains (e.g., Whittaker et al. 2010), because flux and grain size depend on bed shear stress or stream power, which typically increase with uplift rates (Tucker and Slingerland 1996, 1997; Robinson and Slingerland 1998; Dade and Verdeyen 2007). In the same sense, because a warmer climate is typically associated with a higher variability in discharge and larger peak floods (Molnar 2001), we anticipate to see a coarsening of grains because larger floods, through larger stream power, are expected to result in the supply of coarser material (Tucker and Slingerland 1996, 1997; Golden and Springer 2006; Armitage et al. 2011; Chen et al. 2018). However, research on modern lakes in the European Alps also demonstrated that a cooler climate is associated with the supply of coarser-grained material and larger volumes of sediment compared to a warmer climate (Glur et al. 2013). Apparently, the relationships between grain size and climate may not be unique. As a third allogenic driving force, a sediment source consisting of high-grade metamorphic crystalline rocks with a low erodibility (Kühni and Pfiffner 2001a) is expected to produce clasts with a larger size than a headwater region made up of low-grade metamorphic lithologies and/or sedimentary rocks (Kodama 1994; Attal and Lavé 2009). Moreover, the structural fabric of the source bedrock, such as the fracture spacing, was reported to control the initial grain-size distribution in the catchement (Neely and DiBiase 2020). These mechanisms potentially contribute to an allogenic control on supply signals in the sense that grain-size development would be driven by the conditions in the hinterland. Besides such allogenic driving forces, alluvial-fan stratigraphies and grain size might also record the occurrence of autogenic and thus internal mechanisms, as documented by field studies (Mohrig et al. 2000; Harvey et al. 2005; Ventra and Nichols 2014; Hajek and Straub 2017) and numerical models (Parker et al. 1998; Pepin et al. 2010; Mouchené et al. 2017). This is potentially important in the foreland, where a fraction of the supplied material is stored, thereby contributing to the stratigraphic record. Such processes operate stochastically and include, e.g., the dynamic modification of channel widths and energy gradients, which can potentially change the transfer probability of material (Paola et al. 2009; Engelder and Pelletier 2013). Yet field-based observations to date have been particularly restricted to systems evolving in endorheic basins (Ventra and Nichols 2014), to small-scale successions at the channel scale (Dade and Verdeyen 2007) or on a limited number of fan systems (Meek et al. 2020). Therefore, on a basin scale, a documentation of how autogenic signals could be preserved in stratigraphy has not been accomplished yet for basins with known allogenic driving forces.

Besides grain size, another parameter to quantify the sedimentary dynamics of a fluvial system is the relative mobility of grains in transport (Fedele and Paola 2007; D’Arcy et al. 2017; Brooke et al. 2018; Harries et al. 2018). For coarse-grained material where the particle size is larger than 2 mm, Fedele and Paola (2007) proposed a solution to calculate the relative mobility function J, which represents, for a particular grain size, the fraction of material in transport opposed to the material in the substrate. The model is based on the principle that selective deposition of grains with a particular grain size results in the build-up of stratigraphic successions, thereby likely reflecting changes in environmental conditions such as sediment supply and/or water flux (Duller et al. 2010; Brooke et al. 2018). Thus, it allows us to explore how capable the dispersal systems on these alluvial fans were to entrain a grain with a given size through time. Similarly to the dependence of grain size, we expect that the relative mobility of grains may be impacted by both allogenic and autogenic mechanisms.

In the Swiss Molasse Basin (SMB), situated north to the Central Alps (Fig. 1), late Oligocene to early Miocene successions consisting of coarse-grained conglomerate beds are interpreted as deposits of alluvial fans (Stürm 1973; Schlunegger et al. 1996; Kempf et al. 1999). Related deposits have been investigated to reconstruct the transport dynamics at the time the sediments were deposited (e.g., Garefalakis et al. 2024), and to infer shifts in sediment supply rates (e.g., Kuhlemann et al. 2001). In the same sense, changes in the stacking pattern (Allen et al. 1991; Sissingh 1997) of individual sections in combination with shifts in grain sizes have been linked to orogenic processes at the Alpine front and in the Central Alps (Kempf 1998; Schlunegger and Castelltort 2016; Garefalakis and Schlunegger 2018). Yet these studies focused on individual stratigraphic sections recording local patterns in the Alpine tectonics, and an attribution of grain-size signals to changes in the driving forces through time is missing for the scale of the SMB. Therefore, we explore whether the long-term tectonic evolution of the Alps, changes in paleoclimate and shifts in exposed bedrock lithologies left detectable grain-size signals preserved in the alluvial-fan stratigraphies across the SMB. The target sections are located both close to the paleo-apex and in the basin axis, with respect to the Alpine thrust front (Fig. 1), and the sediments cover the time span between c. 31 and 13 Ma.

The aim of this paper is to quantify the sediment dynamics of late Oligocene to early Miocene megafan deposits across the entire SMB that were deposited during a period when allogenic controls on sediment routing changed significantly in the adjacent orogen (Kuhlemann and Kempf 2002; Schlunegger and Kissling 2022). We do so by collecting grain-size data from c. 50,000 grains embedded in conglomerate beds along 15 stratigraphic sections (Fig. 1C), and by calculating the relative mobility of grains in transport. Using the results, we evaluate how sensitive records of alluvial fan stratigraphies are to changes in climate and tectonic conditions over the 18 Myr period. We then discuss whether effects related to allogenic and autogenic mechanisms were preserved in the Molasse conglomerates.

The SMB is made up of a clastic sedimentary wedge that accumulated during Oligo-Miocene times, recording the tectonic and erosional evolution of the adjacent Alpine hinterland (Pfiffner 1986; Allen et al. 1991; Pfiffner et al. 2002; Schlunegger and Kissling 2022). As such, the conglomerate beds and particularly the sizes of grains preserved in these sediments are anticipated to preserve signals that could be related to the development of the Alps. The SMB itself is subdivided into the undeformed and flat-lying Plateau Molasse at distal positions relative to the Alpine orogen, and the tilted, folded, and thrusted Subalpine Molasse adjacent to the Alpine front (Fig. 1B, C). A zone where Molasse thrust sheets dip in the opposite way to each other, referred to as the “Triangle Zone,” marks the transition between the Plateau and the Subalpine Molasse (Fig. 1C).

At the scale of the SMB, the deposits form two large-scale regressive megacycles (Sinclair et al. 1991; Kuhlemann and Kempf 2002), where each cycle comprises a transition from marine to terrestrial conditions (Matter et al. 1980; Pfiffner 1986). The first megacycle records the change, at 30 Ma, from a deep marine basin (North Helvetic Flysch and Lower Marine Molasse) to terrestrial conditions (Lower Freshwater Molasse). Approximately 10 Myr later, the “Burdigalian transgression” marked the start of the second megacycle with deposition of the shallow marine Upper Marine Molasse (Allen et al. 1985; Schlunegger et al. 1997b; Jost et al. 2016; Garefalakis and Schlunegger 2019). Deposition of the Upper Marine Molasse ended with a regression around 16.5 Ma, when terrestrial conditions lead to the formation of the Upper Freshwater Molasse (Kuhlemann and Kempf 2002; Kälin and Kempf 2009). Along the Alpine thrust front, the basin fill comprises several hundreds of meters thick suites made up of conglomerate beds (Fig. 2). They record the occurrence of coarse-grained alluvial megafans, which either interfingered with a fluvial floodplain (Lower and Upper Freshwater Molasse) or a shallow marine environment (Upper Marine Molasse) in more distal positions of the basin (Kuhlemann and Kempf 2002). In general, the conglomerate beds (Fig. 2) are made up of rounded and subrounded clasts embedded in a sandy matrix (Stürm 1973; Bürgisser 1980; Schlunegger and Castelltort 2016; Garefalakis and Schlunegger 2018).

The Central Alps, bordering the SMB to the South (Fig. 1C), evolved in response to the continent–continent collision between the European and African continental plates (Schmid et al. 1996). The collision was accomplished through a rollback subduction of the European plate beneath the African continent (Kissling and Schlunegger 2018). It resulted in the growth of the Alpine topography, starting at c. 32–30 Ma, when the break-off of the European mantle lithosphere slab beneath the Central Alps (Schlunegger and Kissling 2022) resulted in rapid tectonic uplift in the area surrounding the Lepontine dome in the core of the Alps (Figs. 1, 3). This uplift was recorded by an increase in exhumation (Schlunegger and Willett 1999; Boston et al. 2017) and erosion rates in the Alpine hinterland (Sinclair 1997), where the landscape changed from an elevated plateau at c. 30 Ma to a highly dissected terrain after 25 Ma (Schlunegger and Kissling 2015; Garefalakis and Schlunegger 2018). Thereafter, between 25 and 22 Ma, the Central Alps were considered to have reached the highest elevation (Fig. 3; Schlunegger and Norton 2013; Krsnik et al. 2021). Around c. 20 Ma, the uplift and exhumation of the External Crystalline Massifs (Fig. 3), situated north of the Lepontine dome (Fig. 1B, C), was initiated (Herwegh et al. 2017; 2020). At this time, the area surrounding the Lepontine dome (Fig. 1B, C) started to experience the highest rates of exhumation, which was accomplished by tectonic erosion through slip along the Simplon fault zone (Schlunegger 1999; Kuhlemann et al. 2001; Boston et al. 2017). The exposure of lithologies with a higher erosional resistance resulted in a northward shift of the drainage divide and initiated a reorganization of the drainage pattern in the Central Alps (Fig. 3). Thereby, an orogen-parallel drainage network with longer flow paths of the Alpine rivers than before 20 Ma was established (Schlunegger et al. 1998; Kühni and Pfiffner 2001b; Spiegel et al. 2001; Stutenbecker et al. 2019). Ongoing orogenesis during Oligo-Miocene times and the related rollback subduction of the European continental plate beneath the African Plate not only controlled the development of the Alpine orogen as summarized above (Kissling and Schlunegger 2018), it also resulted in the continued deflection of the foreland plate, the formation of accommodation space in the SMB, and the formation of a wedge-shaped basin geometry. In the basin axis, > 10–20 km from the Alpine thrust front, this is reflected by the continuously increasing subsidence rates, on average from < 0.2 mm/yr before 25 Ma to > 0.3 mm/yr after 20 Ma (Schlunegger and Kissling 2015, 2022).

The tectonic development of the Alps, as summarized above, has been considered to have controlled the orogen-scale pattern of sediment supply to the Molasse basin (Hay et al. 1992; Schlunegger 1999; Kuhlemann 2000; Kuhlemann et al. 2001; Schlunegger et al. 2001; Bernard et al. 2021; Schlunegger and Kissling 2022; Garefalakis et al. 2024). Such fluxes, which have been established at the scale of the entire basin, were calculated by previous authors from the preserved mass (Hay et al. 1992) or volume (Sinclair 1997; Schlunegger 1999; Kuhlemann 2000; Kuhlemann et al. 2001; Schlunegger et al. 2001) of sediments in the SMB. Their results show that the slab break-off at 30 Ma and the rapid tectonic uplift of the Lepontine area (Fig. 3) initiated a wave of erosion. This caused an increase in sediment supply from c. 2,000 km³ Myr⁻¹ before 30 Ma, to c. 15,000 km³ Myr⁻¹ thereafter (Hay et al. 1992; Sinclair 1997; Schlunegger 1999; Kuhlemann 2000; Kuhlemann et al. 2001; Schlunegger et al. 2001). Around 25–22 Ma (late Chattian to Aquitanian), when the Central Alps reached the highest elevation (Schlunegger and Norton 2013; Krsnik et al. 2021), the supply of sediment to the SMB occurred at the highest rates of c. 15,000–22,000 km³ Myr⁻¹ (Fig. 3). The post-20 Ma increase in the length of the flow paths in the Alps together with the exhumation of rocks with a high erosional resistance (Bernard et al. 2021), were considered to have caused a decrease of the sediment supply rates to the basin to c. 12,000 km³ Myr⁻¹ (Fig. 3). Together with a eustatic sea level high, this resulted in the Burdigalian transgression (Keller 1989; Schlunegger and Kissling 2024). A re-establishment of terrestrial conditions in the entire SMB at c. 16.5 Ma occurred in response to a drop in the eustatic sea level and a rapid increase in sediment supply to the basin to c. 16,000 km³ Myr⁻¹, after which the supply rates decreased to a constant value of c. 12,000 km³ Myr⁻¹ (Fig. 3).

Garefalakis et al. (2024) provided estimates for long-term sediment fluxes based on two-dimensional fan transects using the self-similar grain-size fining model (Fedele and Paola 2007). For the Rigi–Rossberg dispersal systems (Fig. 1C) approximately 40 km² Myr⁻¹ of material accumulated on the fan between c. 27 and 25 Ma (Fig. 3). The Lake Thun fan system (Fig. 1C) recorded sediment accumulation of c. 16 km² Myr⁻¹ between 25 and 23 Ma (Fig. 3), and c. 6 km² Myr⁻¹ sediments were deposited on the fan of the Hörnli–Töss suite (Fig. 1C), roughly between 15 and 13 Ma (Fig. 3). These supply rates of sediment have been linked to regional tectonic processes that were considered to promote the formation of accommodation space, thereby outpacing a possible influence of a climatic driver (Garefalakis et al. 2024). While this interpretation was made for individual sections, we also explore in this work whether this observation persists at the scale of the entire SMB. We particularly expect that variations in grain size closely follow interpreted trends in sediment supply rates, which were determined by Kuhlemann (2000) and Kuhlemann et al. (2001) using the volumes of preserved sediment (Fig. 3). The reason for this is that both variables should be positively correlated with stream power or shear stress of rivers supplying the material with sources in the orogen (Dade and Verdeyen 2007).

The Central Alps and the SMB experienced various climate changes during their evolution, which could have had the potential to influence alluvial-fan sedimentation in the region. Globally, the period between c. 32 and 27 Ma was characterized by a relatively stable and cold paleoclimate (Fig. 3), as revealed by the relatively high values of the oxygen isotope record (δ¹⁸O; Zachos et al. 2001). This is consistent with paleoclimate interpretations that are based on calcareous seeds of charophytes collected in the SMB (Berger 1992; Schlunegger et al. 2001). Following this period, the global oxygen isotope record indicates a shift towards lower δ¹⁸O values and thus warmer conditions (Zachos et al. 2001), which is referred to as the Late Oligocene Warming Event (LOWE; Fig. 3). Around c. 25.5 Ma, a shift towards a warmer and drier, possibly also stormier (Schlunegger and Norton 2013) climate was interpreted from stable carbon and oxygen isotope data collected in the SMB (Berger 1992; Schlunegger et al. 2001). Following the LOWE, the Oligo-Miocene Transition (OMT; Fig. 3) at c. 23 Ma marks the shift towards a cold period, albeit of short duration only (Zachos et al. 2001). During the mid-Miocene, after a phase of a relatively unstable climate, temperatures increased between c. 17 and 15 Ma. This rise corresponds to the Miocene Climate Optimum (MCO; Fig. 3) and is documented in both the global records (Zachos et al. 2001) and in the carbon and oxygen isotopes collected from pedogenic carbonates in the SMB (Methner et al. 2020; Krsnik et al. 2021). This shift towards a warmer climate in the Molasse basin is also recorded by flora and in related paleoclimate reconstructions (Mosbrugger et al. 2005). Eventually, the MCO was followed by a cooling around c. 14 to 13 Ma (Middle Miocene Climate Transition, MMCT; Fig. 3). This change is recorded globally as well as in the deposits of the SMB (Zachos et al. 2001; Methner et al. 2020; Krsnik et al. 2021; Mosbrugger et al. 2005). Because climate change has the potential to influence the production of sediment through weathering and erosion, it could potentially impact the grain-size distribution of the material supplied to the SMB.

The 15 stratigraphic sections analyzed in this study are situated along the strike direction of the present-day Alpine front (Fig. 1C), and the corresponding deposits span the time interval between 31 and 13 Ma (Fig. 4). They are found in the western, central, and eastern SMB (Fig. 1C). We selected these sections because they have been studied in previous contributions regarding their sedimentological properties and the petrographic composition of the embedded clasts (Büchi 1958; Bürgisser 1981; Matter and Weidmann 1992; Schlunegger et al. 1997c; Kempf et al. 1999; Von Eynatten 2003; van der Boon et al. 2018). Additionally, they have been dated in the framework of numerous paleontological and magnetostratigraphic studies in the western and central (Schlunegger et al. 1996; 1997a; Kälin and Kempf 2009; Engesser and Kälin 2017) and the eastern part of the basin (Kempf et al. 1997; Kempf and Matter 1999; Kälin and Kempf 2009). The original age models are based on correlations of the magnetostratigraphic polarities and biostratigraphic data with the CK95 (Cande and Kent 1995), the Astronomically Tuned Neogene Time Scale ATNTS2004 (Gradstein et al. 2004), or the ATNTS2012 age models (Vandenberghe et al. 2012).

In the west, we analyzed six stratigraphic sections (Figs. 1C, 4). Their thicknesses (Fig. 4) range from 500 m (Emme) to > 3000 m (Lake Thun), and they record the development of the Beichlen–Lake Thun–Napf dispersal system between c. 31 Ma and 14.7 Ma (Fig. 4; Schlunegger et al. 1996). In the central part of the SMB we analyzed two stratigraphic sections (Figs. 1C, 4), which are the c. 1700-m-thick Rigi and c. 1200-m-thick Rossberg suites (both c. 29.5–24.7 Ma; Fig. 4). Here, the Rigi and Rossberg sections comprise the proximal–distal suite of the Rigi dispersal system (Stürm 1973; Garefalakis et al. 2024). In the east, we analyzed seven successions (Figs. 1C, 4). Their thicknesses cover the range between a few hundred meters (Hörnli, Töss) to > 2000 m (Necker; Fig. 4). The sediments encountered in these profiles document the development of the Speer–Kronberg–Hörnli dispersal system between c. 31 Ma and 13 Ma (Fig. 4; Kempf et al. 1999).

The top-most units exposed along the Lake Thun, Prässerebach, Necker, and Rigi sections comprise amalgamated conglomerate beds, which record the occurrence of a paleo-apex or the point-of-entry of the Alpine streams into the foreland basin (e.g., Kempf et al. 1999) (Fig. 4C). Measurements of paleo-flow directions reveal a radial drainage of the corresponding paleo-rivers (Kempf et al. 1999; Schlunegger and Norton 2013; Garefalakis and Schlunegger 2018). In contrast, the deposits encountered at the other sections were deposited near the basin’s axis and thus in a more distal position > 10 km farther downstream of the Alpine thrust front (Kempf et al. 1999; Garefalakis et al. 2024). There, sole marks, imbrications, and cross-beds reveal a NE-directed axial discharge before 18 Ma, and in a NW-oriented dispersal thereafter (Kempf et al. 1999; Kuhlemann and Kempf 2002; Schlunegger and Norton 2013; Garefalakis et al. 2024). These more distal sections largely expose alternations of conglomerate and mudstone beds (Fig. 4C).

Determinations of the petrographic composition of the conglomerates showed that the material deposited between c. 31 and 28 Ma (Emme, Thur, and Steintal sections) are mainly made up of sedimentary clasts (sandstones, siliceous limestones, limestones, and dolomites) with a few crystalline constituents (Schlunegger et al. 1998; Kempf et al. 1999). After that time, the relative abundance of crystalline clasts (granite, gneiss, and rhyolite) started to increase from 10–20% for the deposits between 28 and 26 Ma (Rigi and Rossberg sections; Stürm 1973; Schlunegger 1995) to > 40% for the younger deposits (Matter 1964; Kempf et al. 1999; Schlunegger and Castelltort 2016). These changes in the petrographic composition of the conglomerates (see Appendix B, Fig. B7, for a compilation thereof) have been used to infer the occurrence of a normal unroofing, where successively deeper lithological units of the Alpine edifice with a higher erosional resistance (Kühni and Pfiffner 2001a, 2001b) were exhumed and eroded. Accordingly, if changes in the exposed lithologies in the hinterland alone would be a significant control on grain size, we would expect to see a continuous coarsening of the material in the SMB.

The previously established age models of the analyzed sections were calibrated to time scales that differ in their numerical ages. Therefore, we recalibrated the original local Magneto Polarity Stratigraphies (MPS) to the global Magnetic Polarity Time Scale (MPTS) thereby using the recent Global Time Scale 2020 as a framework (GTS2020; Gradstein et al. 2020). For this, we considered information from micro-mammalian fossil sites (Table B17; Appendix B) and new mapping results in the study area (Hantke et al. 2022; Strasky et al. 2022) to harmonize the age model for all sections. The revised age model for the Lake Thun, Töss, Hörnli, Rigi, and Rossberg sections (including Sattel and Einsiedeln) have been taken from Garefalakis et al. (2024). For each section, the calibration of the magnetostratigraphic polarities to the GTS2020 is given in Appendix B (Figs. B1–B6), and Table B16 indicates which local age model (original or revised MPS) has been considered for the recalibration to the global time scale (GTS2020).

The grain-size distributions (GSD) from the 15 sections were collected on digital photos (Panasonic Lumix FT-5) taken from outcrops of conglomerate beds (see Fig. 4 for number of sites or outcrops per section). We considered nearly all outcrops that were accessible and exposed. The GSD of the Lake Thun, Töss, Hörnli, and Rossberg sections (which includes the data of the Sattel and Einsiedeln sections) were taken from Garefalakis et al. (2024), and the data of the Rigi section was taken from Garefalakis and Schlunegger (2018). On each photograph (3 to 6 per outcrop > 5 m²), taken from an outcrop (or site) at a distance of 1 to 1.5 m, the longest visible axes of 100 grains > 2 mm were measured using ImageJ (Rasband 1997) following a grid-based Wolman point-count approach (Wolman 1954) and the measuring protocol outlined in Garefalakis et al. (2023). Accordingly, we imposed a grid of 10 cm × 10 cm for all photographs (example in Fig. 2F) and measured the longest visible axis underneath a grid point. Because we focus on detecting signals in the coarse-grained fraction of the fluvial deposits, we did not incorporate any sand-GSD. Note that in contrast to Garefalakis et al. (2023) we did not apply an occlusion correction to our grain-size data. As shown by the aforementioned authors, such a correction would lead in a linear shift of the grain-size data only, without changing the outcome of our analysis. Because the sizes of fluvially transported clasts have been shown to follow a logarithmic-, log-normal, or a gamma distribution (Friedman 1962; Church and Kellerhals 1978; Vaz and Fortes 1988; Armitage et al. 2011), we transformed all grain-size data into a natural logarithm-space (ln(GSD)) for our model calculations (cf. Harries et al. 2018).

For each of the 15 sections, we combined all grain size data collected from the various outcrops along a profile (GSD_outcrop) to one bulk database per section referred to as GSD_section. We then calculated the values of the D₅₀ and D₈₄ grain-size percentiles from the GSD_outcrop and the GSD_section databases. Note that, for instance, the D₅₀ percentile values correspond to the grain sizes where 50% of the grains are smaller than or equal to the specific grain size of this percentile. We combined the data into GSD_section because we are interested in gaining information on how grain size and the related transport dynamics evolved at the scale of the entire SMB. In addition, we analyzed the evolution of the GSD_outcrop by calculating moving averages of the grain-size percentiles. For this, we considered a moving age window spanning 25% of the total age range (see, e.g., Fig. 4) and for each window a minimum of two datapoints (i.e., sites or outcrops).

According to the underlying principles of the model, the mobilization of coarse-grained material occurs if grain size-dependent thresholds during bankfull discharge are exceeded (Paola et al. 1992; Bunte et al. 2004; Parker et al. 2007; Wickert and Schildgen 2019). Consequently, for the case that J > 1, grains with a given size are likely in transport, whereas for J < 1, grains with specific sizes are preferentially stored in the substrate and eventually in the stratigraphic record. Thus, in the case where J = 1, a specific grain in a GSD has an equal probability (i.e., 50%) of being either in transport or stored in the substrate (Armitage et al. 2011; Whittaker et al. 2011; Brooke et al. 2018). Consequently, determining the grain size at J = 1 allows the inference of how effectively the dispersal systems on these alluvial fans were able to entrain material with the corresponding size through time.

Here, D_k is an individual grain size, and D¯ and σ are the mean and sample standard deviation of the GSD. We calculated the self-similarity variables from our datasets and tested whether the inference of similarity is statistically valid (see statistical tests below). Because we applied a log-transformation to the GSDs (cf. Armitage et al. 2011; Harries et al. 2018) we refer to the related self-similarity as (ln(GSD)). Similar to the grain-size percentiles, we calculated the ξ₅₀ and ξ₈₄ of the self-similarity distributions for the GSD_outcrop and the GSD_section.

Here, the parameters a_g, b_g, and c_g (all ≥ 0) describe the shape of the relative mobility function for gravel-size particles (subscript g) and characterize the incipient motion of these (Fedele and Paola 2007). In particular, the parameter a_g quantifies the mobility of all grains from a given GSD (if a_g increases, the particle mobility increases in proportion to grain size), b_g influences the rate of the relative mobility (if b_g increases, particles with larger sizes become less mobile than smaller ones), and c_g describes the minimum probability at which grains of all sizes are transported (if c_g > 1, then all grains are likely in transport). In log space, the parameter a_g was shown to lie in the range of 0.30–0.38, b_g between 0.78–1.29 (including uncertainties), and c_g tends to be close to 0 (see also results) but could also be fixed at 0.01 (Harries et al. 2018).

The values for ξ, the mean grain size D¯, and the sample standard deviation σ are based on the GSD_section and are calculated through bootstrapping (BT, section below), thereby using the log-transformed data.

We used standard statistical tests including: i) the Shapiro-Wilk (SW) test for normality (Shapiro and Wilk 1965; Shapiro et al. 1968) to check if our grain-size datasets follow a normal distribution, and ii) the Kolmogorov-Smirnov two-sample (KS2) test (Hodges 1958) to explore if the GSDs are similar to each other. Both tests were performed at a significance level of α = 5%, equivalent to the 95% confidence level. We note that the SW test is accurate if the sample size does not exceed 5000 (Shapiro and Wilk 1965). Because four sections exceed this number (i.e., Rossberg, Prässerebach, Lake Thun, and Necker; see results), we additionally report the SW statistics for these (Appendix A). In these cases, the data likely follow a normal distribution if the SW statistics tend to approach 1. To perform the SW test, we simulated a normal distribution based on the average and standard deviation of the ln(GSD_section) datasets for each section and tested each outcrop datasets (ln(GSD_outcrop)) against these. For the sake of completeness, both tests (K2S and SW) were performed on the values of the log-transformed (ln(GSD)) datasets, as well as on the self-similarity variables (ξ(GSD) and ξ(ln(GSD)).

For all calculations including grain-size data, such as the relative mobility function, we applied a bootstrapping approach. For this, we calculated the 95% C.I. upon resampling 100 grain-size values with replacement to simulate 10⁴ scenarios. We applied this procedure to both the GSD_outcrop and GSD_section datasets, and to the related ξ values. The application of BT also accounts for potential biases from sample number and statistical outliers in the GSDs.

We measured a total of 49,900 grains along the 15 stratigraphic sections. The number of measured grains per section ranges between 1000 and 7300 (Fig. 5; i.e., 100 grains per outcrop). The D₅₀ and D₈₄ values of the outcrops (GSD_outcrop) range from c. 20 to 105 mm and from 34 to 190 mm, respectively. When plotted against the depositional ages of the sediments, these D₅₀ and D₈₄ values (Fig. 5) display a large scatter within related uncertainties (i.e., 95% C.I. of individual percentiles). The results of the moving average show a narrower spread, with the D₅₀ ranging from c. 27 to 91 mm and the D₈₄ values between c. 48 and 175 mm. At the scale of individual sections, most of the analyzed successions do not reveal clear grain-size trends through time, as neither a consistent fining nor a coarsening is observed or could be accommodated by the spread of the uncertainties (Fig. 5). This is particularly the case for those conglomerates that were deposited in the axial drainage. Yet we note that for the Prässerebach, Lake Thun, and Necker sections, which record deposition at the Alpine thrust front close to the paleo-apex, a shift towards larger grains occurred between c. 25 and 22 Ma (Fig. 5). When the data are merged into the GSD_section, the same grain-size percentiles of interest range from c. 30 to 50 mm for the D₅₀, with an overall long-term average of 40 mm, whereas the values of the D₈₄ range from c. 60 to 100 mm, scattering around an average of 80 mm (Fig. 6). Importantly, except for the deposits around 24 Ma, these values are nearly constant through time and disclose neither a clear temporal nor a spatial trend. Any variations are within related uncertainties (Fig. 6).

The percentiles (ξ₅₀ and ξ₈₄) of the self-similarity variables (GSD_section) are nearly identical irrespective of the depositional age and location of the sediments (Fig. 7). As such, the ξ₅₀ ranges approximately between –0.01 and 0.1, and the ξ₈₄ varies roughly between 0.9 and 1.2 (Fig. 7). In addition, both percentile values are within the associated uncertainties (95% C.I.) of the (GSD_section). Therefore, for all sections and through time, the values that characterize the variations in grain sizes and thus the mean and the standard deviations of the underlying GSD_section are proportional to each other. This suggests that also the sorting properties and underlying transport processes of the dispersal systems on the individual fans were likely similar. Moreover, this remarkable similarity persists despite the different ages of the explored deposits, and the different positions of the sections and of the sampling sites in the basin (Fig. 1C). It also persists despite the different petrographic composition of the conglomerate beds (e.g., the relative abundance of crystalline versus sedimentary clasts) and thus despite the inference that the source lithology of the individual dispersal systems changed through time (Fig. B7 in Appendix B).

The cumulative frequency curves of the individual section’s GSDs display a large spread between each other, and they all reveal a right-skewed shape (Fig. 8A). The logarithmically transformed grain-size data disclose evenly distributed cumulative GSDs for both, the upper and lower tails (ln(GSD)), Fig. 8B). This transformation reveals that the skewed tails are not a result of the sampling approach, but that the GSD datasets follow a log-normal distribution (cf. Harries et al. 2018). The cumulative frequency curves of the self-similarity variables ξ(GSD) (dimensionless; Eq. 2) all collapse on the same curve, yet the left-skewed distributions still disclose that the underlying GSDs indeed follow a log-normal distribution (Fig. 8C). The calculations of the self-similarity variable using the log-transformed grain-size data, i.e., (ln(GSD)), returns cumulative frequency curves that match the simulated normal distribution, and the skews in the tails are nearly absent (Fig. 8D). The behavior of these cumulative frequency curves is observed not only for the data collected from individual outcrops (GSD_outcrop), but also for the merged dataset representing an entire section (GSD_section) and for the corresponding transformations (see Figs. A3–A7 in Appendix A).

We applied the KS2 test to explore whether the grain-size data are statistically similar to each other across the individual sites or outcrops for each section, and we applied the SW test to check if the log-transformed grain size and self-similarity variable datasets follow a normal distribution. Both tests were performed at the 95% confidence level (Figs. A8–A11 in Appendix A for the KS2 and SW). For both, the GSD_outcrop and ln(GSD_outcrop) datasets, the statistical similarity is on average c. 51% (KS2 in Fig. 8A, B), while the transformation of the GSD_outcrop into the self-similarity variable ξ reveals that both, the ξ(GSD_outcrop) and ξ(ln(GSD_outcrop)) datasets are on average c. 99% similar across the individual sites of each section (KS2 in Fig. 8C, D). The SW test result of the ln(GSD_outcrop) and the ξ(ln(GSD_outcrop)) datasets are on average c. 75% (SW in Fig. 8C, D). For sections where > 5000 grain sizes were measured, the SW statistics are on average > 0.97 for the log-transformed datasets (ln(GSD_outcrop)) and ξ(ln(GSD_outcrop)), and thus close to 1 in all cases (Appendix A). Transforming the GSD into the self-similarity variable ξ largely improves the statistical similarity between these datasets (KS2 test results; Fig. 8C, D), while a logarithmic transformation of both datasets markedly improves the probability that the data are normally distributed (SW test results; Fig. 8B, D). Therefore, for the calculations of the relative mobility function, which requires the data to follow similar distributions (similarity) and implicitly assumes a normal distribution of the data, we considered the (ln(GSD_outcrop)) datasets only (Fig. 8D). Please note that we used the bulk dataset, i.e., (ln(GSD_section)), of each stratigraphic section for the corresponding calculations.

The median best-fit solutions of f(ξ)_{(med. BT)} (solid curves in Fig. 9; Appendix A, Eq. A4) are plotted against the ξ_{(med. BT)} (dotted curves in Fig. 9; Eq. 2) calculated for the ξ(ln(GSD_section)) distributions as probability density curves. The modeling results yield a high similarity between the model outcomes across all sections, as each of the individual best-fit solutions collapse on the same curve (Fig. 9; see Appendix A, Fig. A12, for results of each section). In more detail, for all sections, the peaks of the f(ξ)_{(med. BT)} (solid) curves of our modeled outcomes are situated slightly above the peak frequencies of the ξ_{(med. BT)} (dotted curves; Fig. 9). Likewise, the left tails and right tails of the median curves of our modeled f(ξ)_{(med. BT)} outcomes are also situated above those of the ξ_{(med. BT)} curves (Fig. 9). In contrast, at the positions of the left shoulders and the right shoulders, the f(ξ)_{(med. BT)} curves of our model outcomes tend to yield smaller frequencies than the corresponding values based on the ξ_{(med. BT)} curves (Fig. 9). Overall, our modeled f(ξ)_{(med. BT)} curves are in good agreement with the spread of the ξ(ln(GSD_section)) distributions (Fig. 9). Please also note that these best-fit solutions depend on the selection of the model input variables, for which we used values as proposed in the literature, and related calculations based on these (C₁, C₂, and C_V, see Figs. A1 and A2 in Appendix A). For instance, for the values for C₁, we used uniformly distributed values between 0.60 and 0.90 (e.g., Fedele and Paola 2007).

The modeling of the best-fit solutions yields in a specific combination of the relative mobility parameters a_g, b_g, and c_g, which characterize the incipient motion of gravel-size particles. Our results fall between 0.24 and 0.27 for a_g, 0.68 and 0.78 for b_g, and are close to 0 (i.e., 2 × 10⁻¹⁶) for c_g (median values; Fig. 10). When plotted against increasing average ages of the deposits exposed along the target sections, the parameters a_g and b_g do not vary within the given uncertainties (95% C.I.; Fig. 10) and do not disclose clear trends. A comparison of our results with those of Harries et al. (2018; also calculated in logarithmic space) revealed that a_g and b_g are c. 15 to 30 and 25% smaller, respectively, thereby revealing a lower probability of larger grains to be entrained. Similarly to their study, the parameter c_g is close to 0, thereby revealing a very limited probability of the coarser grain-size fraction to be transferred in transport. Lastly, the model parameters’ numerical values depend (amongst other variables; see above) on the spread of the underlying GSD, which are coarser-grained in Harries et al. (2018) compared to our study.

The shape of the relative mobility curve (Eq. 3) is largely controlled by the numerical values of the model-input variables (C₁, C₂, and C_V, see Appendix A, Fig. A2) and the three parameters, a_g, b_g, and c_g (see above and Appendix A, Fig. A3). As such, higher values of the three parameters in general tend to decrease the relative mobility of particles with increasing grain size. Based on the model outcomes for f(ξ) and their fitting to the ξ(ln(GSD_section)) data, which yielded the results for the related parameters, a_g, b_g, and c_g, we calculated the values for the relative mobility function J (Eq. 3) as a function of the ξ(ln(GSD_section)) distributions for all sections. Our results reveal a congruent pattern for the curves of J(ξ(ln(GSD_section))) (Fig. 11A). The uncertainties expressed by the 95% C.I. largely overlap with the curves characterizing the individual sections. This shows that the solutions for the relative mobility function, when displayed as curves and for all sections, collapse on the same distribution (Fig. 11A), whereas the results of previous studies displayed some variations between different systems, yet within related uncertainties only (Brooke et al. 2018; Harries et al. 2018). In particular, for the case where J = 1, the median values of the (ln(GSD_section))-distributions range between –1.97 and –1.80 with an overall average of c. –1.89 (inset in Fig. 11A and Figs. A16 and A18a in Appendix A).

The concept of relative mobility offers a possibility to explore how environmental conditions are reflected by the temporal evolution of the grain-size data. For this, the non-dimensional (ln(GSD_section)) curves (Fig. 11A) need to be back-transformed into the dimensional grain size values (Eq. 4). Our results reveal a pattern where the median curves generally overlap within the 95% C.I. for the various sections (Fig. 11B). For the intersection where J = 1 (i.e., the “equally mobile” grain size where a grain has the probability of 50% to be in transport or stored in the substrate) the corresponding median critical grain-size values vary between c. 8 and 16 mm with an overall average of 11 mm (inset in Fig. 11B and Figs. A17 and A18B in Appendix A). The cross-comparison of these results revealed that all median critical grain-size values show some variations through time, but scatter within the 95% C.I. (Fig. 12), thereby revealing a congruent pattern for the various sections through time. In addition, the intersection of the relative mobility curves (Fig. 11B) with the median grain-size values of the D₅₀ and D₈₄ (e.g., Fig. 6) yield direct information about the relative mobility of these specific grain-size percentiles. The results revealed that these percentiles have on average a c. 75% (for the D₅₀) and 88% (for the D₈₄) larger probability of being stored in the substrate (i.e., J < 1) than being entrained (i.e., J > 1; Fig. 13A, B). Therefore, the relative mobility of gravel-size particles at all sections is relatively low but, considering the spread of the uncertainties, broadly constant through geological time.

Our results show that the rivers feeding the fans preferentially entrained grain sizes smaller than 11 mm on average (Fig. 12), whereas coarser-grained material had a high probability of being stored in the substrate. In particular, the values of the critical grain sizes for J = 1, which contain information on the sedimentary dynamics related to selective entrainment and deposition, are markedly smaller than the corresponding section’s D₅₀ and D₈₄ values (Fig. 6). Therefore, the D₅₀ and D₈₄ have higher probabilities to be stored in the substrate (75% and 88%, respectively; Fig. 13) than being entrained, indicating a relatively low mobility of the gravel in these rivers. We interpret that the relatively low probability of entrainment of the D₅₀ and D₈₄ is due to the role of intermittent high-discharge events controlling the entrainment of material on alluvial fans during relatively short periods (Dury 1961; Molnar 2001; Navratil et al. 2006; Leenman et al. 2022). Such bankfull or high discharge events are more likely to control the motion of the coarser grain sizes such as the D₈₄ or the D₉₆ (Ferguson 2007; Philips and Jerolmack 2019). In our case, we calculated a transport probability of c. 10% for a D₈₄ size clast (Fig. 13). This is in agreement with the results of Begin and Inbar (1984), who estimated a 10% probability for the movement of particles with sizes related to the D₉₀. This is also consistent with studies that have calculated the transport probability on recent alluvial fans where sediment transport has been dominated by intermittent floods and sheet flows. There, the D₈₄ is c. 50–60% (fans in Death Valley, USA; Brooke et al. 2018) or the mean clast size is c. 60–70% (fans in the Iglesia Basin, Argentinian Andes; Harries et al. 2018) more likely to be stored in the substrate than in transport.

If we focus on the suite of amalgamated conglomerates that were deposited close to the paleo-apex and, thus, at relatively proximal sites on the corresponding fans, we find evidence for a supply control on grain-size signals. In particular, sections which start at their bottom with distally deposited conglomerate and mudstone alternations, and end with a proximal suite of amalgamated conglomerate beds towards their top (e.g., Prässerebach, Lake Thun, Necker; Fig. 4C), show trends of increasing grain sizes through time (Fig. 5). These coarsening trends largely persist if a moving average is calculated for the data (Fig. 5). They are particularly observed for c. 25 and 22 Ma old conglomerates, falling into a time when sediment fluxes to the basin were increasing to the highest rates and when the relative abundance of crystalline clasts increased in the conglomerate beds (Fig. 3). This increase in sediment supply rates was interpreted as a response to Alpine tectonic events either at the thrust front (e.g., Pfiffner 1986; Kempf et al. 1999) or in the core of the Alps (Schlunegger and Castelltort 2016). It also occurred at the same time when: i) the crystalline rocks were exhumed to the surface (Schlunegger et al. 2001), ii) the Alps reached their highest elevations (Schlunegger and Norton 2013; Krsnik et al. 2021) and, when iii) global climate experienced a transient phase towards cooler conditions (OMT, Fig. 3). We thus envisage an interpretation where the coarsening observed at these sections reflects the combined effect of a tectonic, topographic, and environmental driving force (Fig. 14). Note that these signals are better visible for the larger-grained percentiles (e.g., D₈₄ values) than for smaller ones (Fig. 14). We interpret this observation using the results of the relative-mobility calculations. They imply that the Molasse stratigraphy preferentially records the occurrence of high intermittent floods with the consequence that the coarse-grained tail of a GSD is better correlated to environmental parameters than, e.g., the D₅₀ and finer-grained percentiles (Whittaker et al. 2011; Schlunegger and Norton 2013; Brooke et al. 2018; Garefalakis and Schlunegger 2018; Schlunegger et al. 2020). We thus observe evidence suggesting that a fan apex has the potential to record a grain-size response to allogenic driving forces in a stratigraphic succession. Additionally, the shift towards amalgamations of conglomerate beds towards the tops of these sections suggests that the sedimentary facies also evolved in response to the aforementioned driving mechanisms. This suggests that facies might be at least as sensitive to allogenic driving forces as grain size alone (Weltje et al. 1998).

For the conglomerates that were deposited in the axial drainage > 10 km farther downstream of the apices, we observe no particular coarsening trends up-section, or, if they exist, they appear to be random (Fig. 5). This is additionally corroborated by the GSD_section grain-size values, which scatter within related uncertainties without disclosing clear trends (Fig. 6). Interestingly and more important, sections recording such facies show no correlation to sediment supply rates, changes in paleoclimate, Alpine tectonic processes (Fig. 14), as well as changes in the relative abundance of crystalline clasts in the conglomerates, and the slightly increasing subsidence rates through time. Therefore, it appears that the downstream distance from the Alpine thrust front was large enough for autogenic mechanisms to buffer or dilute any (allogenic) supply-controlled grain-size signal from the source area (which appear to be recorded close to the apex; see paragraph above), thereby averaging grain sizes over stratigraphic timescales of 10⁵ to 10⁶ years. Autogenic mechanisms have the potential to result in modifications of energy gradients and channel widths, thereby changing the transport probability of grain sizes throughout the system (Engelder and Pelletier 2013). We note that autogenic processes were considered as occurring over short time spans only (i.e., < 10⁴ years; Einsele et al. 1991), yet recent studies have shown that these self-organized mechanisms can also act over longer time scales (> 10⁶ years; Hinderer 2012; Hajek and Straub 2017). Indeed, if autogenic processes operated continually at short timescales, as is likely, they may shred signals over the long term, with the observed long-term record potentially representing an extended history of short-term autogenic dynamics (Jerolmack and Paola 2010; Wang et al. 2011; Paola 2016; Romans et al. 2016; Li et al. 2018; Straub et al. 2020). Moreover, the lag time of an environmental signal to arrive in the foreland basin could further amplify these effects (Castelltort and Van Den Driessche 2003; Schlunegger and Castelltort 2016; Li et al. 2018). Accordingly, in more distal sites of a foreland basin, signals related to allogenic forces might be preserved in the stratigraphy only if their amplitude is much larger than that of autogenic processes, or if the response time in the basin to tectonic processes in the source area is relatively short (Romans et al. 2016; Li et al. 2018; Sharman et al. 2019).

We analyzed long-term records of grain-size trends offered by fluvial conglomerate beds in the Swiss Molasse Basin (SMB) through time and the relative mobility of clasts. For this, we focused on 15 sections recording alluvial megafan sedimentation in the SMB during c. 18 Myr from Oligocene to Miocene times (c. 31 to 13 Ma). At these sections, we measured the sizes of c. 50,000 grains, and we statistically analyzed the resulting grain-size distributions. We applied the concept of self-similarity and relative mobility to this dataset, and we estimated the critical grain size of a particle that is equally likely in transport or preferentially stored in the substrate (J = 1). Because tectonic processes in the source area, the supply rates of the sediment to the basin, the climatic conditions at the time the material was deposited, and changes in the lithology and thus erodibility of exposed bedrock in the Alps changed significantly during the time when gravel sedimentation occurred, we anticipated marked differences in the grain size and relative clast mobility through time.

Our study reveals that at the scale of the SMB, the grain-size values, such as the D₅₀ and D₈₄, were nearly constant through time, despite the substantial changes in the study area mentioned above. If changes in supplied grain size did occur, they were too small to be statistically significant within the stratigraphic record. This is particularly the case for conglomerates that were deposited in the basin axis > 10 km downstream of the thrust front. Likewise, the relative mobility of grains on the dispersal systems did not significantly change through time, showing that the sorting processes and mechanisms related to selective deposition and transport were similar through geological time. This is also valid for the relative mobility of grains at the scale of the individual sections, where potential differences are accommodated by the spread in the uncertainties. However, for those sections where the most proximal facies of a fan is preserved, particularly the coarser grain-size percentiles (e.g., D₈₄ in our study) do record shifts in the environmental conditions, as previously proposed in literature.

In conclusion, the nearly identical sizes of grains preserved in the fluvial conglomerate beds deposited in the axis of the Swiss Molasse Basin (SMB) points to the effect of autogenic self-organized mechanisms on a fan, which might shred possible grain-size-related supply signals and prevent them from being recorded in the stratigraphy. In addition, the nearly constant relative mobility of grains, which is also inferred from the grain-size distributions in the conglomerate beds, implies that the mechanisms of selective transport and deposition were almost identical during the 18 Myr of record we document in the SMB. The data also point towards a low relative mobility of grains through Oligocene and Miocene times when the conglomerates were deposited. Yet, we find that allogenic signals, such as shifts in supply rates, might be preserved for those systems that recorded the most proximal facies (such as the apex of a fan), particularly if the larger percentiles and the scale of individual sections are considered. Therefore, the proximal facies of fluvial conglomerates, often preserved in the uppermost part of stratigraphic sections, likely are the best places to observe records of grain-size signals and transport signals attributed to allogenic processes, whereas autogenic mechanisms dilute or buffer these signals at the more downstream facies as illustrated in this dataset.

The study has been conceptualized by PG, FS, and ACW. PG carried out the fieldwork, collected the samples, and analyzed the data with additional scientific input by AHdP, ACW, DM, and FS. PG prepared the manuscript and the figures with contributions from FS, ACW, DM, and AHdP.

Philippos Garefalakis reports financial support provided by the Swiss National Science Foundation (SNSF) (Grant No. P1BEP2_200189). Ariel Do Prado’s salary was funded by the Marie Sklodowska-Curie Innovative Training Network S2S (Grant No. 860383).

The authors would like to thank editor Dustin E. Sweet and associate editor Doug Edmonds for handling and John Southard for editing our manuscript. We further would like to thank Melissa Lester for managing the manuscript. We are also grateful to Laure Guerit and Alexander B. Neely for their detailed and constructive reviews, which greatly improved this work.