Temporal variations in rockfall and rockwall retreat rates in a deglaciated valley over the last 11 ka

To calculate the power law predicted retreat rate, power law fit parameters must be defined. Previous studies have shown that above a certain magnitude, the frequency of rockfalls and landslides appears to decay as a power-law (Strunden et al., 2015; Barlow et al., 2012; Dussauge et al., 2003; Guzzetti et al., 2003; Hunger et al., 1999). Below this magnitude, which is commonly referred to as “rollover”, there is no clear trend. The rollover might be due to undersampling of small rockfalls by data collection methods (Lim et al., 2010; Stark and Hovius, 2001) or may have a physical explanation reflecting different slope failure processes (Stark and Guzzetti, 2009; Guthrie and Evans, 2004). To fit a power law distribution of the form:


Power-law Predicted Rockwall Retreat Rate
To calculate the power law predicted retreat rate, power law fit parameters must be defined. Previous studies have shown that above a certain magnitude, the frequency of rockfalls and landslides appears to decay as a power-law (Strunden et al., 2015;Barlow et al., 2012;Dussauge et al., 2003;Guzzetti et al., 2003;Hunger et al., 1999). Below this magnitude, which is commonly referred to as "rollover", there is no clear trend. The rollover might be due to undersampling of small rockfalls by data collection methods (Lim et al., 2010;Stark and Hovius, 2001) or may have a physical explanation reflecting different slope failure processes (Stark and Guzzetti, 2009;Guthrie and Evans, 2004). To fit a power law distribution of the form: to the observed data, first the rollover value V o has to be determined before the scaling exponent ƅ can be obtained. Events smaller than the rollover value V o are excluded in the power law fit. To estimate V o , two approaches are used: in the first approach, we compute nonlinear fits to the data for all possible rollover volumes. From these, the best fit as determined by the highest R 2 value is selected (blue line in Fig. DR1). The second approach follows the procedure described by Clauset et al. (2007). In this approach, the choice of the rollover volume and its corresponding maximum likelihood estimate (MLE) for the scaling exponent ƅ are tested using the Kolmogorov-Smirnov (KS) statistic, in which the lowest value describes the best fit (orange line in Fig. DR1). Since both approaches are reasonable but result in quite different rollover volumes (0.29 vs. 1.07 m 3 , respectively), further investigation is necessary. Figure DR1 shows two peaks in R 2 , corresponding to two levels in the KS value. This implies deviations from pure power law 2 behavior in the observed frequency-magnitude dataset and creates a conflict between intentions to find the best fit results and to retain data. To illustrate, at a rollover volume of 0.29 m 3 , 30% of recorded rockfalls are discarded in contrast to 64% at a rollover volume of 1.07 m 3 . We further investigate this issue by studying the power law fit results over this range of rollover volumes. Since the position of the lowest KS value is governed mostly by data noise, we substitute the upper limit by the position of onset of the second level in KS and R 2 peak, 0.92 m 3 . Results shown in Figure DR2 reveal a systematic gradient in the scaling exponent b and derived rockfall projections with higher rollover volumes.
For the computation of the scaling exponent ƅ, we compare nonlinear least squares regression fits with MLE, as suggested by Clauset et al. (2009) (Fig. DR2, panel a). We avoided using the linear least squares regression methods in log space as they are known to be sensitive to deviations in the tail of the distribution and to incorrectly underestimate the uncertainty associated with the data (Catani et al., 2016;Clauset et al., 2009). As Figure DR2 shows, the nonlinear regression's 95% confidence interval is substantially smaller than the ML uncertainty calculated after Clauset et al. (2007) and almost fully contained within the latter's confidence band. Furthermore, only the ML uncertainty captures the spread inherent in the data, as illustrated by the systematic gradient depending on the rollover volume. Consequently, we consider the nonlinear regression results as too constrained and report result values in the rest of this study as interval of ML estimates for the aforementioned range of rollover values.
The total long-term eroded volume was obtained by integrating the frequency density of rockfall magnitudes multiplied by the event volume following Barlow et al. (2012). The integration boundaries were chosen based on the minimum rollover volume (V min = 0.29 m 3 ) and the estimated volume of largest historic event ~130 years ago (V max = 10 4 m 3 ). Based on the total eroded volume, the corresponding intermediate timescale (power law predicted) wall retreat rate was calculated.

Long Timescale (11 kyr) Rockwall Retreat Rates
Talus deposits have been used to measure Holocene rockwall retreat rates (Hinchliffe and Ballantyne, 1999;Rapp, 1960;Hoffmann and Schrott, 2002;Siewert et al., 2012). The magnitude of rockwall retreat is calculated by measuring the talus volume, compensating for the density difference between talus bulk and intact bedrock, and dividing by the rockwall surface area. To calculate rockwall retreat rate, the rockwall retreat magnitude is divided by the talus production time (i.e., timespan elapsed since deglaciation). In this study, we selected five talus deposits that have well-constrained lateral extents and are not incised by waterfalls (Fig. 1). However, only three taluses (WT2, WT3, and ET1) were considered in the calculation of Holocene rockwall retreat rate. Talus WT1 and WT4 yield biased retreat rates due to complex morphology and fluvial incision, respectively.

3
The talus volume is estimated by reconstructing the lateral extent and surface of the talus slopes with TLS data in the JRC software and using shaded relief derived from 2-meter resolution digital elevation model swissALTI3D. Point cloud data corresponding to talus slopes were compiled and filtered to remove vegetation using the curvature filter. This procedure produced a smooth talus surface. Using the topographic meshing tool and an ordinary kriging algorithm, triangle meshes of smoothed talus surfaces were constructed and used in volume calculations. Since the Lauterbrunnen Valley contains up to 550 meters of sediments (Kellerhals and Haefeli, 1985), but only the upper 35 meters have been deposited since the Younger Dryas (Bodmer et al., 1973), the present-day height of the valley floor was set as a boundary between the present-day talus slope and the underlying valley fill when calculating talus volumes. Fig.  DR3 summarizes the steps involved in TLS data processing for quantifying talus volumes. Point cloud data are archived in opentopography.org (see: https://doi.org/10.5069/G9GQ6VX7).
We used a mean density correcting factor of 0.77 to compensate for the density difference between talus and intact rock. This value assumes a talus bulk debris density of 2.0 g/cm 3 and an intact rock density of 2.6 g/cm 3 based on data collected from settings with similar lithologies (e.g., Krautblatter et al., 2012 andSass andWollny, 2001).
The talus accumulation timespan or the time that has elapsed since deglaciation in the Lauterbrunnen Valley is assumed to be at most 11 kyr. This is based on radiocarbon ages of moraines located in the upper Lauterbrunnen Valley, ~5 km south of Stechelberg (see Fig. 1). These ages indicate that this area was ice-free at 10.39 ± 0.15 ka (Wipf, 2001). Additionally, 10 Be surface exposure ages of boulder and bedrock samples from a neighboring trough valley, Hasli, imply that deglaciation was completed in this area at the end of the Younger Dryas at 12.2-10.8 ka (Wirsig et al., 2016).
Rockwall retreat rates calculated from talus volumes are only rough estimates due to the aforementioned assumptions. The largest source of uncertainty is the unknown bedrock profile beneath each talus and the transition from talus deposit to valley infill. To constrain these, and to reduce possible errors, talus deposits with complex morphology (e.g., WT1), and those that have been significantly modified (e.g., WT4) are removed from wall retreat calculations. Additionally, we have calculated variations in the long-term rockwall retreat rates for different scenarios in which talus volumes and production time are altered systematically and the corresponding rockwall retreat rates are compared (Table DR4). Increasing talus volumes by 10%, 20% and 50% yield rockwall retreat rates of 0.36, 0.40, and 0.50 mm/yr, respectively. Similarly, using talus age of 2 ka, 5 ka, and 7 ka, we calculate rockwall retreat rates of 0.52, 0.73 and 2 mm/yr, respectively. Therefore, in this study we report the average long-term retreat rate obtained from talus volumes to be > 0.33 mm/yr. This value is based on maximum talus age and minimum talus volume.

Extrapolated Return Period for Large Events
Assuming that the observed volume distribution follows a power law, the distribution can be extrapolated outside the volume range and time window (Dussauge et al., 2002). Following the approach used in Strunden et al. (2015), MLE results are extrapolated. The occurrence of a rockfall event of 10 3 m 3 is predicted to be every 3.34-6.40 years. These return periods, however, do not mean that an event larger than 1000 m 3 will reliably occur every 3-6 years, as it is only an average over long periods. Though these predictions agree with the previously estimated return periods by Strunden et al. (2015), there is a large uncertainty associated with them. This is due to the large extrapolation of rockfall volumes from our dataset (largest event 267.27 ± 4.39 m 3 ) to an event of 10 3 m 3 . Furthermore, the rockfall distribution law may vary over time. The observation duration of 5.2 years is not much longer than the event recurrent period, and therefore, the temporal variability of rockfall distribution cannot be fully assessed. In addition, rockfall activity may depend on loading conditions such as climate, leading to temporal fluctuations of rockfall activity and power law exponents. A larger time window for rockfall activity may be achieved by utilizing historical datasets, provided they are accurate and complete.