Articles | Volume 8, issue 3
Earth Surf. Dynam., 8, 579–593, 2020
Earth Surf. Dynam., 8, 579–593, 2020

Research article 08 Jul 2020

Research article | 08 Jul 2020

The impact of earthquakes on orogen-scale exhumation

The impact of earthquakes on orogen-scale exhumation
Oliver R. Francis1,2, Tristram C. Hales1,2, Daniel E. J. Hobley1, Xuanmei Fan3, Alexander J. Horton1, Gianvito Scaringi4, and Runqiu Huang3 Oliver R. Francis et al.
  • 1School of Earth and Ocean Sciences, Cardiff university, Main Building, Park Place, Cardiff, CF10 3AT, UK
  • 2Sustainable Places Research Institute, Cardiff University, 33 Park Place, Cardiff, CF10 3BA, UK
  • 3The State Key Laboratory of Geohazard Prevention and Geoenvironment Protection, Chengdu University of Technology, Chengdu, 610059, Sichuan, China
  • 4Institute of Hydrogeology, Engineering Geology and Applied Geophysics, Faculty of Science, Charles University, Albertov 6, 128 46 Prague 2, Prague, Czech Republic

Correspondence: Oliver R. Francis (


Individual, large thrusting earthquakes can cause hundreds to thousands of years of exhumation in a geologically instantaneous moment through landslide generation. The bedrock landslides generated are important weathering agents through the conversion of bedrock into mobile regolith. Despite this, orogen-scale records of surface uplift and exhumation, whether sedimentary or geochemical, contain little to no evidence of individual large earthquakes. We examine how earthquakes and landslides influence exhumation and surface uplift rates with a zero-dimensional numerical model, supported by observations from the 2008 Mw 7.9 Wenchuan earthquake. We also simulate the concentration of cosmogenic radionuclides within the model domain, so we can examine the timescales over which earthquake-driven changes in exhumation can be measured. Our model uses empirically constrained relationships between seismic energy release, weathering, and landsliding volumes to show that large earthquakes generate the most surface uplift, despite causing lowering of the bedrock surface. Our model suggests that when earthquakes are the dominant rock uplift process in an orogen, rapid surface uplift can occur when regolith, which limits bedrock weathering, is preserved on the mountain range. After a large earthquake, there is a lowering in concentrations of 10Be in regolith leaving the orogen, but the concentrations return to the long-term average within 103 years. The timescale of the seismically induced cosmogenic nuclide concentration signal is shorter than the averaging time of most thermochronometers (>103 years). However, our model suggests that the short-term stochastic feedbacks between weathering and exhumation produce measurable increases in cosmogenically measured exhumation rates which can be linked to earthquakes.

1 Introduction

Surface uplift of a mountain range is controlled by the balance of additive uplift processes and removal of material by surface, typically fluvial, processes. Earthquakes produce rock uplift and equally importantly generate exhumation via landsliding (Avouac, 2007; Keefer, 1994; Marc et al., 2016a). Landsliding events can scour deeply into bedrock, causing many hundreds or thousands of years of exhumation of the bedrock surface in a geologically instantaneous moment (Fig. 1) (Li et al., 2014; Marc et al., 2019; Parker et al., 2011; Stolle et al., 2017). Existing mass balances on single (Parker et al., 2011) or sequences (Li et al., 2014, 2019; Marc et al., 2016b) of earthquakes demonstrate that landslide volumes of large thrust earthquakes are comparable to, and may exceed, rock uplift. However, earthquakes are not the only rock uplift process in mountain belts; aseismic mechanisms such as viscous and elastic crustal deformation (Meade, 2010; Simpson, 2015), lithospheric delamination (Hales et al., 2005; Molnar et al., 1993), and isostatic rock uplift (Molnar, 2012; Molnar et al., 2015) can also contribute. The total time-averaged mass balance, and any net surface uplift, is also affected by erosion between earthquakes (Hovius et al., 2011; Marc et al., 2019; Yanites et al., 2010). Therefore, the contribution of earthquakes to the generation of long-term surface uplift of mountains remains poorly constrained (England and Molnar, 1990; Li et al., 2019).

Figure 1Landsliding can generate large volumes of sediment but only the largest consistently transport sediment off the hillslope. The volume of bedrock eroded is also dependent on the thickness of regolith on the hillslope before the landslide occurs. Panel (A) shows the 2016 Aranayake landslide in Sri Lanka which was triggered by intense rainfall in heavily weathered soils. The failure occurred along the regolith–bedrock interface and very little bedrock was eroded. Panel (B) shows that in bedrock-dominated areas, such as the Usoi Dam in Tajikistan (produced by the 1911 earthquake), significant volumes of regolith can be generated by landsliding but little measurable erosion as the regolith produced has remained in the catchment as a landslide dam for over 100 years. The 2008 Wenchuan earthquake produced over 60 000 landslides (C), many of which have significantly smaller transport lengths than the hillslope (highlighted in red) and so cannot be easily removed from the catchment. For these reasons, we define landsliding as weathering rather than as erosional agents. Panel (A) is sourced from (last access: 3 July 2020), (B) is from (last access: 3 July 2020), (C) is sourced from Google Earth, using imagery provided by Landsat/Copernicus and Maxar Technologies, 2019 and draped over a digital terrain model. © Google Earth.

Despite the importance of earthquake-triggered landslides in the total erosion budget of mountain belts, there is little evidence of large earthquakes in most sedimentary or exhumation records. Increased rates of sedimentation and changes in sedimentary characteristics, linked to the abundance of loose sediment, have been identified in lakes and reservoirs immediately proximal to large faults (Howarth et al., 2012; Zhang et al., 2019). However, these pulses in sedimentation can be difficult to separate from extreme hydrological events (Zhang et al., 2019). Steep slopes within mountain belts, as well as typical observations of plentiful exposed bedrock on those slopes, have been used to support the idea that earthquake-generated sediment is rapidly removed from orogenic belts (Dingle et al., 2018; Li et al., 2014; Niemi et al., 2005; Parker et al., 2011). However, observations collected after earthquakes suggests that regolith, transportable sediment produced by landsliding or weathering, remains in low-order catchments (Pearce and Watson, 1986) (Fig. 1a), single large landslide deposits (Korup and Clague, 2009; Stolle et al., 2017), landslide dams (Fan et al., 2012; Ouimet et al., 2007) (Fig. 1b), and valley fills (Blöthe and Korup, 2013) for up to 104 years. An estimated 44 % of the sediment produced in the Himalaya is stored in the long term in some form before it exits the mountain range, demonstrating the importance of storage and sediment recycling to orogen-scale sediment fluxes (Blöthe and Korup, 2013). Similar findings are found in studies on post-earthquake sediment fluxes: in New Zealand, up to 75 % of the sediment produced by the 1929 Murchison earthquake remained in catchments 50 years after the earthquake (Pearce and Watson, 1986), while in Taiwan between 92 % and 99 % of the sediment produced by the Chi-Chi earthquake remained in catchments 8 years after the earthquake (Hovius et al., 2011). The slow removal time of earthquake-derived sediment suggests residence times could exceed recurrence times of earthquakes in tectonically active mountain ranges. If significant volumes of sediment remain in the landscape for times longer than the recurrence period of earthquakes, large earthquakes could contribute significantly more to surface uplift than currently assumed (Li et al., 2014; Marc et al., 2016b; Parker et al., 2011). Constraining the long-term contribution of earthquakes to a mountainous orogen is impossible without fully constraining the rate at which coseismic landslide deposits are removed from the mountain range and the contribution of non-seismic erosion and uplift. The rock uplift minus the coseismic landsliding produced by earthquakes is commonly assumed to be equal to surface uplift, which directly assumes there is no sediment cover in the mountain range, and thus all landslides are in bedrock (England and Molnar, 1990; Li et al., 2014, 2019; Marc et al., 2016b). If sediment remains in orogens for extended periods of time, it is likely it will be remobilised, reducing the erosion and increasing the surface uplift of an earthquake. Remobilisation of coseismic sediment can also occur between earthquakes, altering sediment fluxes and estimates of erosion rates from cosmogenic nuclides (Andermann et al., 2012; Dingle et al., 2018; Yanites et al., 2009). Constraining any aseismic uplift and how it affects coseismic uplift is also be important for understanding whether earthquakes build mountain ranges (Hovius et al., 2011; Li et al., 2014; Parker et al., 2011; Royden et al., 1997). Observations of the 2008 Wenchuan earthquake and its aftermath give us insight into the export of landslide-generated regolith from mountain belts. Over 60 000 landslides were produced by the earthquake, with a total volume of approximately 3 km3 (Li et al., 2014) (Fig. 1c). In mountain ranges, there are significantly more small landslide deposits than large ones, such that the magnitude and frequency of these volumes follows a power law (Malamud et al., 2004a; Marc et al., 2016a; Stark et al., 2001). Only the largest, or most mobile, landslides are most likely to deposit sediment directly into the channel; for the Wenchuan earthquake, less than half of the total volume of the regolith produced is connected to the channel network immediately after the earthquake (Li et al., 2016) (Fig. 1c). Mapping of the landslide deposits reveals much of the landslide material produced by the earthquake remains on the hillslope and in small-order channels, and many deposits are undistributed by erosional processes (Fan et al., 2018; Zhang et al., 2016). When looking at an orogen as a whole, landsliding only transports sediment a very small distance, and hence we define landsliding as a weathering rather than erosional process.

The rate at which landslide deposits can be evacuated from a catchment is typically related to the capacity of the fluvial system to transport the influx of sediment (Croissant et al., 2017; Yanites et al., 2010). If a significant proportion of sediment produced by an earthquake is not able to be eroded by the fluvial system, the residence time of the sediment is likely to be increased due to the need for stochastic hillslope processes, such as debris flows, to remobilise the sediment into the channel before it can be exported (Bennett et al., 2014; Croissant et al., 2019; Hovius et al., 2011). The timing of debris flow triggering is related to the interaction of storms and slope hydrology, which cannot be easily predicted, while the volume remobilised by debris flows is primarily controlled by the non-linear process of sediment entrainment during run out (Horton et al., 2019; Iverson, 2000). The rapid stabilisation of landslides after earthquake, combined with the stochastic triggering of debris flows, could be one reason why large pulses of sediment associated with single earthquakes are rare or absent from downstream sinks (Fan et al., 2018; Zhang et al., 2019). Rather than a sink receiving a large impulse of sediment, which can be easily recorded via a change in average grain size or sedimentation rate, the rate change is instead smeared or shredded by stochastic sediment transport processes like flooding or debris flows (Jerolmack and Paola, 2010). The averaging time for different measures of erosion rate (e.g. cosmogenic vs. thermochronometric) may strongly affect the probability of measuring a single earthquake. If the recording time of erosion (103 years for cosmogenic radio nuclides or 104−5 years for thermochronometry methods) is similar to or less than the return period of large earthquakes, then any difference between short- and long-term erosion rates could be due to the influence of earthquakes (Kirchner et al., 2001; Ouimet, 2010). By investigating the variation of erosion rates with varying timescales, or with coseismic landslide density (Niemi et al., 2005; West et al., 2014), we may be able to identify the long-term impact of earthquakes on orogen-scale erosion rates.

In this paper, we use a zero-dimensional volume balance model to explicitly track earthquake generated sediment through time in a hypothetical orogen based upon the Longmen Shan region. We use the tracked sediment thickness in order to understand how earthquake-triggered landslides (EQTLs) affect exhumation and surface uplift at orogenic scales. Our model co-varies the amount of aseismic uplift in the orogen, imposed earthquake magnitude–frequency relationships, and both the timing and maximum magnitude of earthquakes, under multiple possible uplift regimes, in order to fully investigate the role of earthquakes in orogen-scale volume budgets. We then use these scenarios to investigate whether earthquakes can be identified from erosion or exhumation records over different timescales and by modelling cosmogenic radionuclide concentrations of sediment leaving the orogen. Finally, we test our hypothesis that the variation in recorded erosion rates in a landscape can be an indicator of seismic activity using a global database of cosmogenic radionuclide derived erosion rates.

Definitions of terms

This paper is precise in its use of terminology around the changes in elevation of the various surfaces we discuss. These follow standard modern definitions (England and Molnar, 1990), but as we are explicitly measuring the generation of regolith and the movement of two surfaces (the topographic surface and the bedrock surface), the potential for ambiguity requires us to clearly define these terms:

  • Erosion is the transport away of material from a site and thus a change in the topography. In our context, it describes the full evacuation of regolith material from the model domain, i.e. removal of sediment from the entire mountain range.

  • Weathering is the in situ conversion of rock into regolith. In our model, rock must be converted to regolith before it can be transported. We explicitly separate the role of EQTLs generating regolith by weathering bedrock from their role as (inefficient) eroders of bedrock and regolith – i.e. we separate their role in producing loose material from their role in transporting that material. EQTLs on average occur on hillslopes near ridge tops, typically with transport lengths less than hillslope lengths, with only the largest impinging on the fluvial system (Li et al., 2016).

  • Regolith is the mobile transportable layer of sediment at the surface. In the model, regolith can be created by two distinct weathering mechanisms: landslides cutting into bedrock to create transportable debris, and soil production by physiochemical processes.

  • Uplift is the increase in elevation of a material or surface in an absolute frame of reference. We distinguish rock uplift, topographic surface uplift, and bedrock surface uplift. Rock uplift is the expected increase in the bedrock surface before considering erosion, it is either produced by coseismic or aseismic means. Bedrock surface uplift is the vertical distance moved by the bedrock surface after erosion. Topographic surface uplift is the total vertical distance moved after considering changes in the bedrock surface and the thickness of the regolith layer.

  • Exhumation is the approach of the rock mass towards the topographic surface and/or the bedrock surface, in the frame of reference of that surface. Our model tracks both surfaces; therefore, it is possible to have a rock uplift event that causes exhumation relative to the bedrock surface without exhumation relative to the topographic surface, by thickening the regolith.

  • Denudation is used almost as a synonym for exhumation, but where the frame of reference is the rock mass or the bedrock surface, and the topographic surface moves towards it.

2 Methods

2.1 Zero-dimensional volume valance model

Here, we present a generalised zero-dimensional mountain volume balance model which we use to test the impact of regolith storage on bedrock surface uplift and exhumation. In the absence of sufficient empirical evidence on the long-term spatial distributions of rock uplift, exhumation, and regolith volumes in mountain ranges, we simulate these parameters by treating the evolution of a landscape as a series of dimensionless seismic volume balances.

In our model, we define the change in topographic surface elevation (ST) with time as

(1) d S T d t = U - E ,

where U (units of length/time) is the thickness of rock entering the orogen during a time step and resulting in rock uplift (calculated as the volume entering the orogen across the area of the model per unit time), while E is the rate of regolith removed (length/time) from the topographic surface and thus is the long-term erosion rate of the orogen. Rock can enter the orogen in two ways: either via shortening and thickening of the crust during coseismic thrusting earthquakes (Uco) or through one of a number of aseismic uplift mechanisms (Uas), such as lower crustal flow (Royden et al., 1997) or lithospheric delamination (Hales et al., 2005). Coseismic deformation (Bonilla et al., 1984; Wells and Coppersmith, 1994), and hence the volume of rock uplifted, scales as a function of earthquake magnitude (Li et al., 2014; Marc et al., 2016b). The addition of mass to a column of crust by thickening will produce an isostatic compensation which will reduce the overall surface uplift response to rock uplift. In our model, we apply a simple compensation based upon the relative densities of the crust and mantle to account for the isostatic response (Densmore et al., 2012; Molnar, 2012; Turcotte and Schubert, 2002). The calculated response is applied immediately to the volume balance and the surface uplift. We do not consider interseismic strain as an uplifting mechanism in this model due its limited contribution to permanent surface deformation (Avouac, 2007). There is ongoing debate to how much surface uplift can be attributed to aseismic vs. coseismic sources and how they interact (Hubbard and Shaw, 2009; Royden et al., 1997). Acknowledging the complexity of the debate, we simplify the aseismic component of uplift and generalise it as the proportion of uplift that cannot be accounted for by Uco. Hence, topographic surface uplift can be represented as

(2) d S T d t = U co + U as - E ,

where the ratio between coseismic uplift rate (Uco) and aseismic uplift rate (Uas) is defined by the term α which represents the proportion of the total uplift rate that is caused by aseismic uplift, such that (1-α)U=Uco and αU=Uas.

In our zero-dimensional model, the thickness of regolith (R) removed from the surface of the orogen is defined as

(3) d R d t = CLRP d t + W d t - E ,

where CLRP (coseismic landslide regolith production) is the average thickness (length/time) of weathered material generated by coseismic landslides, all of which is assumed to be transportable, and W (length/time) is the thickness of rock weathering caused by all other mechanisms (simply the thickness of material removed from the bedrock when there is no regolith cover). W is included to ensure that erosion can continue even when regolith is not present. In our model, weathering does not occur when there is a covering of regolith as the background weathering rates for our study site in the Longmen Shan region are unknown. This way of including weathering in our model allows it to be an emerging property rather than a fixed rate. The rate of elevation change (length/time) of the bedrock surface (SB) can now be described as

(4) d S B d t = U co + U as - CLRP d t - W d t .

We can now define our rate of surface uplift again as a combination of the bedrock surface elevation and the regolith thickness:

(5) d S T d t = d S B d t + d R d t .

The model represents the average topographic surface uplift, regolith generation, and bedrock surface lowering for the area (A) of coseismic displacement of the largest possible earthquake for a fault found within a mountain belt. The length of the modelled area (L) is set by the length of the surface rupture on the fault that generates the maximum earthquake, while the width is the distance to the estimated line of zero strain based upon the dip of the modelled fault (θ) and the focus depth (D) (Li et al., 2014).

(6) A = L D θ

As surface uplift rates for mountain ranges are hard to define (England and Molnar, 1990), we set the model to a flux steady state, where U is set to equal the long-term erosion rate (E). For each time step in the model, an earthquake with Mw>5 is randomly chosen from a power law distribution, and the coseismic rock uplift volume associated with this earthquake is calculated using an empirical scaling relationship (Li et al., 2014) between magnitude and rock uplift volume:

(7) log V u = 1.06 ( ± 0.22 ) M w - 8.40 ( ± 1.44 ) ,

where Vu is the volume of rock uplift generated by an earthquake of magnitude Mw. This volume is scaled by α and divided by the model area A to calculate Uco. Mw 5 earthquakes are the smallest that regularly produce coseismic landsliding so represent the smallest earthquakes of interest to our study (Marc et al., 2016b). We use an optimising algorithm to fit the uplift produced by Eq. (7) to ensure the model remains in a flux steady state. The algorithm uses the uncertainty within Eq. (7) to fit the model parameters so that the average uplift produced during a time step is equal to the long-term erosion rate. The use of the aseismic uplift scaling parameter α has the effect of increasing the time-averaged rock uplift of time steps of small earthquakes and decreasing the rock uplift of large earthquakes.

Regolith is generated in the model based on calculations of bedrock lowering by CLRP and weathering by other mechanisms. We use the scaling between landslide volume (Vl) and earthquake magnitude proposed by Malamud et al. (2004b) (Fig. 2a)

(8) log V l = 1.42 M w - 11.26 ( ± 0.52 ) ,

Vl is converted to a depth of landsliding by dividing by the area of the model space (A). The area of the landscape affected by landsliding of the largest earthquakes is greater than the model space, so a scaling is applied based on the landslide density of the Wenchuan earthquake. Alternative models of coseismic landslide volume as a function of seismic moment (Marc et al., 2016a; Robinson et al., 2016) cannot easily be scaled into a zero-dimensional model space due to their reliance upon earthquake source depth and landscape metrics. These models describe the relationship between earthquake magnitude and total landslide volume as a curve around a hinge magnitude. The shaking produced by an earthquake correlates with the length and width of its surface rupture; however, the width (depth) of the rupture is limited by the thickness of elastic crust. At a maximum magnitude (Mw 6.75), the scaling between earthquake magnitude and shaking changes resulting in a curved relationship between total landslide volume and earthquake magnitude. This has the effect of reducing the importance of large earthquakes in the surface uplift balance. As our chosen model does not include this threshold, the larger earthquakes of this model will be more erosive than other models. All empirical models relating coseismic landslide volume and earthquake magnitude have large uncertainties in them. We acknowledge these uncertainties by applying a normal distribution of error using the uncertainty bounds on the landsliding volume produced by each earthquake (Fig. 2a), reducing the difference between the different models. The total new regolith generated by coseismic landslides is then calculated as the difference between the average depth of landsliding and the average thickness of the regolith.

(9) CLRP = V l A - R

Figure 2The average depth of regolith produced by an earthquake is impacted by the earthquake magnitude and the thickness of regolith that is on the hillslope before the earthquake occurs. (a) Scaling of landslide volume with magnitude from Malamud et al. (2004b) and the average residence time of landslide sediment in the Longmen Shan region based upon an erosion rate of 0.62 mm yr−1. The two lines represent the minimum and maximum volumes of landsliding generated, within the bounds on Eq. (8). (b) Variability of regolith production, expressed as volume per area, with existing depth of regolith on the hillslope, for four representative earthquake magnitudes. Coloured areas represent the variability of the landslide volume produced by an earthquake, randomly sampling within the bounds of Eq. (8). These coseismic landslide regolith production functions (CLRPFs) emerge from the model rather than being set in advance, and the variability at each magnitude is driven by noise inherent in the relationship between magnitude and landslide size (Eq. 8). The inset shows the probability distribution function for regolith thickness across the whole model run, integrating the effects of the CLRPF through time. A small but non-zero spatially averaged modal regolith thickness persists, but significantly larger thicknesses regularly occur. (c) Variability of surface elevation through time for model runs with identical earthquake sequences but with varying additional proportions of aseismic uplift. The assumption of steady state prevents any long-term permanent uplift, and so all variations in surface elevation are driven by the sequence of earthquakes and changes in regolith thickness. Increasing the aseismic contribution to uplift reduces the uplift of large earthquakes, resulting in much less variable surface uplift and therefore exhumation.


2.2 Model implementation: Longmen Shan

We test our model using the Longmen Shan region due to the wealth of studies of the area both prior to and after the 2008 Wenchuan earthquake. These studies allow for the small number of parameters in our model to be well constrained. The Longmen Shan region marks the eastern margin of the Tibetan Plateau and the western edge of the Sichuan Basin (G. Li et al., 2017). Hillslopes are at their threshold steepness with pervasive bedrock and limited channel storage (Li et al., 2014; Ouimet et al., 2009; Parker et al., 2011). The front of the range is dissected by three parallel dextral-thrust oblique-slip thrust faults, two of which, the Yingxiu-Beichuan and Pengguan faults, ruptured during the 2008 Wenchuan earthquake (Densmore et al., 2010; Liu-Zeng et al., 2009). Prior to this earthquake, geodetic measurements recorded limited shortening rates suggesting a possible role for lower crustal flow in driving surface uplift (Clark et al., 2005; Kirby et al., 2003; Royden et al., 1997). However, significant shortening associated with the Wenchuan earthquake supports an important, possibly exclusively coseismic, surface uplift element (Hubbard and Shaw, 2009).

To apply the model to the Longmen Shan region, we use a power-law relationship of earthquake frequency and magnitude derived from historical earthquake records (Z. Li et al., 2017):

(10) N > M w = 3.93 - 0.91 M w .

N is the number of earthquakes that occurs above a certain magnitude in a year. The smallest earthquake we model is a Mw 5, which occurs on average every 5 years. This relationship gives a return time of 1816 years for earthquakes of the same magnitude as the Wenchuan earthquake. Other studies (Densmore et al., 2007; G. Li et al., 2014, 2017) have proposed a return time of anywhere between 500 and 4000 years for a Mw 7.9 earthquake. We use the uncertainty in the frequency of Wenchuan-sized earthquakes to vary Eq. (9) to test the impact of earthquake frequency on exhumation and surface uplift. We use an apatite fission track-derived exhumation rate of 0.62(+0.14-0.08) mm yr−1 (G. Li et al., 2017) to represent the long-term sediment flux from the orogen (E). The model area is set by Eq. (6) using parameters derived from the Wenchuan earthquake. The length (L) is the surface rupture of the Wenchuan earthquake (240 km), the focal depth (D) was between 14 and 18 km deep, and the dip angle ranged between 40 and 90 (Li et al., 2014). These parameters yield an area (A) of 2600 km2. We run the model for 25 Myr to allow multiple analyses over various timescales and vary α between 0 and 1.

2.3 Exhumation calculations within the model

We calculate exhumation for 2 kyr intervals, which is the average time our model takes to exhume 1.2 m of rock through the rock–regolith interface given our assumed value of E. This depth is chosen as it broadly representing the recording timescale for cosmogenic radionuclides (Gosse and Philips, 2001). Exhumation is calculated as the difference between rock uplift (Uco+Uas) and bedrock surface uplift (SB) over the recording time. We randomly choose 10 000 2 kyr samples from each 25 Myr run to produce a distribution of exhumation rates. We investigate the change in exhumation rates due to different proportions of coseismic and aseismic rock uplift and varying earthquake frequency and maximum earthquake magnitudes.

2.4 Cosmogenic radionuclide calculations

We also calculate the cosmogenic 10Be flux out of the model through time. For each time step, we add 10Be atoms to the system using published production rates and attenuation lengths to simulate the depth profile of cosmogenic concentrations (Balco et al., 2008; Braucher et al., 2011; Granger and Muzikar, 2001).

(11) P ( z ) = P 0 e - z ( ρ / Λ )

The production rate (P) of 10Be decreases exponentially with depth below the topographic surface, z, from the production rate at the topographic surface (P0) based upon the density of the bedrock (ρ) and the attenuation length (Λ). The production rates (atoms g−1 yr−1) and attenuation length (g cm−2) depend on the radiation being modelled. In our model, we simulate spallation (production rate of 5.784 atoms g−1 yr−1, attenuation length of 160 g cm−2), fast muons (production rate of 0.0418 atoms g−1 yr−1, attenuation length of 4320 g cm−2), and slow muons (production rate of 0.014 atoms g−1 yr−1, attenuation length of 1500 g cm−2), and combine them to give a total concentration at depth intervals set by the long-term erosion rate (Braucher et al., 2011; Granger and Muzikar, 2001). When an earthquake generates regolith, the top depth intervals are mixed, and the constant erosion rate is applied to remove regolith from the surface. After a spin-up time of 10 kyr, the model tracks concentration of 10Be in the sediment leaving the model. The spin-up time is the time required for the concentrations to reach a stable value which is perturbed by earthquakes. As erosion in the model is constant, and set to the long-term exhumation rate, any change in concentration represents the effect of stochastic EQTLs on exhumation.

3 Results and discussion

3.1 Coseismic landslide regolith production

Within our model, regolith generated by the largest earthquakes can remain on hillslopes for ∼1000 years (Fig. 2a). The average thickness of new regolith that is produced in an earthquake (expressed as volume per area; i.e. a depth) is a strong function of both the pre-existing depth of regolith prior to new failures, the magnitude of the earthquake, and stochastic differences in the volume of landslides for a given earthquake magnitude (Fig. 2b). The primary control on the total landslide volume produced by an earthquake magnitude is the strength of the shaking, with topography and rock strength as secondary factors (Marc et al., 2016a; Valagussa et al., 2019). If shaking can produce similar volumes of landsliding regardless of how much bedrock or regolith is on the hillslopes, landslide deposits in a mountain range with widespread regolith will contain less fresh bedrock, as regolith will make up a greater proportion of material mobilised by the earthquake. As regolith makes up a greater volume of the landslide deposits, less bedrock weathering occurs. This effect will be particularly powerful in areas where large earthquakes occur frequently in similar locations. However, the distribution of earthquake magnitudes exerts a stronger influence on total regolith production through time, as more frequent small earthquakes can only ever weather small depths of regolith from the bedrock (Fig. 2b). The decline of coseismic landslide regolith production (CLRP) with existing regolith thickness is reminiscent of soil production functions described for soil mantled landscapes (Heimsath et al., 1997), and by analogy, we term the non-linear relationship between regolith production rate and the average depth of weathering by landslides seen here a “coseismic landslide regolith production function” (CLRPF). However, unlike a “traditional” soil production function, two elements of stochasticity are inherent to a CLRPF. One reflects the role of shielding of the bedrock surface (SB) from lowering when the regolith layer is thicker than the average depth of the generated landslides (Larsen et al., 2010), and that thickness is dependent on the past history of landsliding in the model. The other reflects the inherent randomness in the size and distribution of the landslides that occur in response to an earthquake of a given magnitude, i.e. within Eq. (8).

As expected, total seismic regolith production is dominated by the largest earthquakes, which produce the largest mean landslide volumes (Malamud et al., 2004b) (Fig. 2b). Summing through time, earthquakes produce 42 % of the total regolith generated by the model; Mw>7 earthquakes account for ∼65 % of the total earthquake-generated regolith and thus 27 % of the total regolith production. However, because smaller earthquakes occur often but produce little regolith, allowing the layer to thin, the time-averaged and spatially averaged regolith layer is predominantly thin – the modal thickness is only 0.02 m, and the mean is 0.03 m (Fig. 2b, inset). Thus, the model shows that although mountain belt regolith cover appears thin almost all the time, at times it can be thick enough to severely affect the short-term exhumation rates of the mountain range. The Longmen Shan region is primarily classified as a bedrock landscape with little storage, but significant volumes of sediment remain in the mountain range after the earthquake, in a similar way to that simulated by the model (Fan et al., 2018; Zhang et al., 2016). Variability in surface uplift through time (Fig. 2c) is affected by the pre-earthquake regolith thickness and therefore the sequence of earthquakes which occur before it. Where large earthquakes are closely spaced in time, pre-existing regolith can limit weathering of the bedrock surface, encouraging uplift of the topographic surface. In catchments close to active faults, the bedrock is likely to be heavily fractured and the shaking is more intense, producing larger landslides with greater densities (Marc et al., 2016a; Meunier et al., 2013; Valagussa et al., 2019). If the regolith is not fully removed from these catchments in between earthquakes, it is possible that the CLRPF may encourage greater surface uplift. In our model, the regolith produced by earthquakes is spread evenly across the landscape, which does not occur in reality. Even in the most impacted catchments in Wenchuan, landslide density is rarely above 10 % km−2, suggesting that remobilisation of landslide regolith on the hillslope may be rare unless the regolith can remain in the catchment for multiple earthquake cycles (Dai et al., 2011; Marc et al., 2015; Parker et al., 2011, 2015). The remobilisation of previous coseismic landslide regolith is, likely to be a local effect mainly impacting catchments close to active fault belts. Ultimately, this interaction between surface uplift and regolith depth is controlled by (1) the time between earthquakes, (2) the magnitude of the previous earthquake, and (3) the rate of regolith removal. The closer together, in both time and space, large earthquakes occur and the slower regolith is removed from hillslopes, the greater the impact of the CLRPF on the surface uplift of a mountain range.

3.2 Regolith generation and volume budgets of earthquakes

Our model demonstrates that the contribution of earthquakes to the uplift and weathering budgets of mountains varies with earthquake magnitude, frequency, and the relative contribution of aseismic weathering and erosion, i.e. erosion or weathering not directly related to earthquakes such as rainfall-triggered landsliding or chemical weathering. When coseismic rock uplift is the dominant uplift mechanism, we can classify four distinct landscape response styles at different earthquake magnitudes (Fig. 3a and c, zones 1–4). On average, an earthquake of magnitude Mw<5.6 lowers the topographic surface (ST). Here, the erosion of regolith (E) out of the model space is greater than the rock uplift produced by the earthquakes (zone 1 in Fig. 3a). For earthquakes with magnitudes 5.6<Mw<7.6, the bedrock surface (SB) rises because the rock uplift rate is greater than the typical regolith generation rate (zones 2 and 3). The regolith thins because the change in the regolith thickness due to CLRP is less than the erosion rate (E) out of the model (zone 2). Conversely, the change in regolith thickness due to CLRP exceeds erosion rate in earthquakes with magnitudes Mw>6.4, so the regolith layer increases in thickness (zones 3 and 4). The largest earthquakes with Mw>7.6 lower the bedrock surface due to the large volumes of regolith produced (zone 4). However, much of the regolith is not removed before the next earthquake, resulting in a net topographic surface uplift, primarily due to thickening of the regolith. For a theoretical purely aseismic uplift scenario, where earthquakes produce landslides but do not create rock uplift (earthquakes produce only horizontal motion, this scenario is not realistic but acts as an extreme end member), earthquakes with Mw>6.5 would cause bedrock surface lowering, while smaller earthquakes would permit bedrock surface uplift due to low CLRP. Earthquakes with Mw>6.5 produce a thick layer of regolith which can persist until the next earthquake, limiting bedrock surface weathering and resulting in net uplift of the bedrock surface.

Figure 3Interplay of changes to the modelled rock uplift rate, the topographic and bedrock surface uplift, and the resulting regolith thickness through time, classified according to earthquake magnitude. The total bedrock and topographic uplift produced by each earthquake magnitude through the model run is summed up and divided by the run time to produce a rate. Panel (a) represents a run with 100 % coseismic uplift, while panel (b) is purely aseismic. Each time a recorded surface intersects the horizontal axis, we separate the chart into a zone which is further described in the text and in panel (c).


Figure 4Kernel density plots (bandwidth of 0.01) of exhumation and denudation rates in various scenarios. Dashed lines indicate the position of the mean ± the standard deviation of the distribution. Each curve is made up of 10 000 samples taken from their respective model runs. (a) Exhumation rates in different uplift regimes. (b) Denudation rates while varying the maximum earthquake magnitude in a run; the run with a maximum magnitude 5 has only earthquakes of a magnitude 5. (c) Denudation rates while varying the frequency of Wenchuan size earthquakes from every 500 to every 4000 years.


Figure 5Variability of (a) topographic surface uplift and (b) the recorded concentration of cosmogenic nuclides leaving the orogen after a representative magnitude 8 earthquake within the model run. Red lines are the mean for the elapsed time since the earthquake, while the grey lines are the real-time concentrations.


3.3 Earthquakes and exhumation

We explore how earthquakes affect exhumation through direct calculation of exhumation of rock at the bedrock surface (SB) relative to the topographic surface. There is very little (∼1 %) variability in erosion rates in model runs with maximum earthquake magnitudes of Mw<5.0. However, the introduction of stochastic weathering of many tens of cm of the bedrock surface by earthquakes with Mw>7 introduces variability in exhumation. When large earthquakes are present, exhumation rates have a standard deviation of 0.055–0.081 mm yr−1 (9 %–14 % of the long-term exhumation rate) and a range of 0.77–0.89 mm yr−1, with the lower figures reflecting a greater contribution of aseismic uplift (Fig. 4a and b). Decreasing the return times of Wenchuan-like earthquakes from 4000 to 500 years produces more variable distributions of exhumation rates, with the standard deviation of exhumation rate increasing from 0.044 mm yr−1 (7 % of long-term average) to 0.076 mm yr−1 (12 % of the long-term average) (Fig. 4c). Taken together, our model results suggest that up to 14 % of the variability in a sample of exhumation rates from a single geographical region could be associated with the time since the last large (Mw>7.0) earthquake. However, this variation may only be seen in exhumation or surface uplift records with recording times of less than 1000 years (Fig. 5a). Pre-Wenchuan earthquake measurements of cosmogenically derived erosion rates are between 40 % and 60 % lower than low-temperature thermochronometrically derived exhumation rates (Ouimet, 2010). Stochastic exhumation of low-concentration bedrock by EQTLs may explain some of that difference.

Cosmogenic radionuclides provide a record of potential earthquake-driven changes in exhumation because they have a relatively short averaging time that is close to the frequency of large earthquakes in many mountain belts. Our modelling results demonstrate the scale of stochastic variability in surface uplift and exhumation. We extend this analysis by simulating cosmogenic concentration in the model to estimate the potential impact of a large earthquake on both the cosmogenic concentration through time and the distribution of cosmogenic concentrations that are likely to be measured. We assume each earthquake thoroughly mixes the regolith down to its average weathering depth. The cosmogenic analysis (Fig. 5b) shows that immediately after a large earthquake, mixing of low-concentration bedrock material with higher concentration regolith lowers the concentration of radionuclides exiting the model. Regolith exiting the mountain range has a lower cosmogenic nuclide concentration for 200–300 years after the earthquake; after this period of low concentrations, there is a peak of concentration higher than the long-term average before a rapid return to the long-term average concentration (Fig. 5b). In the case of a representative magnitude Mw∼8 earthquake, the concentration falls initially by around one-third. However, the process of mixing also increases the concentration of nuclides close to the regolith–bedrock interface compared to the values before mixing, so that as the regolith is slowly eroded through time the lower half releases concentrations greater than the long-term average. Therefore, in landscapes with frequent Mw>7 earthquakes and regular long-term storage of regolith, it is possible to record more variable cosmogenic concentrations than might be expected, including positive as well as negative excursions from the long-term mean (Fig. 5b).

Figure 6Reanalysis of detrital cosmogenic radionuclide derived denudation rates for mountain belts around the world compiled by Harel et al. (2016), in the context of peak ground acceleration and tectonic environment. (a) A box plot indicating the median (central orange line), quartiles (end of box) and the range (“the whiskers”) of denudation rates in the analysed localities ordered by their average seismicity (in brackets), defined as the maximum peak ground acceleration (PGA) of a 475-year return period. Points indicate values outside the range, ±1.5 times the interquartile range. (b) Standard deviation (SD) of denudation rates for each mountain belt against seismicity represented by the 475-year return PGA (m s−2). (c) Standard deviation of denudation rates for each mountain belt compared to the standard deviation of steepness indexes. The areas of highest variability are found in steep, tectonically active mountain ranges.


These modelling results provide testable predictions of the exhumation-related changes to cosmogenic concentration caused by large earthquakes. Hence, we look to examine whether the predicted variability might represent some of the variability associated with cosmogenic erosion rates in seismically active areas using a global dataset compiled by Harel et al. (2016). Harel et al. (2016) collated detrital cosmogenic 10Be concentrations for 59 geographical areas, separated into areas of similar climates, and recalculated the erosion rates using consistent production rate and shielding corrections. We limit our sampling to those sites with >18 measurements and basins larger than 105 m2 to limit sampling bias (Dingle et al., 2018; Niemi et al., 2005). We compare the probability density distribution of erosion rates from within those geographic regions to seismicity, as represented by the 475-year return peak ground acceleration (PGA) derived from a global seismic hazard map (Giardini et al., 1999; Harel et al., 2016) (Fig. 6a). Due to the size of the geographical areas, there may be multiple seismic hazard levels recorded; we simply use the mean value to classify the area. The use of a single number to characterise a large area can underestimate the potential PGA. While a single number may not accurately describe the entire area, it does allow us to compare the variability of denudation rates with seismic hazard. We crudely classify the regions as dominantly coseismic or dominantly aseismic: regions with thrust faults (identified from a literature review) and erosion rates greater than the median are deemed coseismic, while slowly eroding regions with no thrust faults are aseismic. “Mixed” regions are those that do not fall under either of these classifications. Coseismic landscapes, as expected, have higher average cosmogenically determined erosion rates with higher standard deviations (Kruskal–Wallis H-test statistic of 14.14; p value of 0.00017) (Fig. 6a and b). Relief is a major control on erosion rates, with steeper catchments having higher erosion rates than shallower ones (Montgomery and Brandon, 2002). The most seismically active mountain ranges are also among the most varied in relief as they have some of the steepest catchments in the world. Therefore, we need to test whether variability in erosion rates is more closely related to variation in catchment steepness or the seismicity of the area. Within the database compiled by Harel et al. (2016), they include a normalised channel steepness index which we use to compare the impact of catchment steepness on erosion rate. The channel steepness index equation is defined by

(12) M x = E K A 0 m 1 n ,

where Mx is the steepness index, E is the erosion rate, K is the erodibility coefficient, A0 is a reference area of 1 m2, and m and n are empirical constants. The index is normalised by assuming fixed values for m and n (Harel et al., 2016). We would expect that in areas with highly variable steepness indexes the erosion rates are also more variable. We find that while areas with higher seismicity have more variable erosion rates, the variation in erosion rates correlates much more closely with the variation in steepness index (Fig. 6b and c). Steeper basins in tectonically active mountain ranges are more susceptible to coseismic landsliding (G. Li et al., 2017; Marc et al., 2016a) and thus will have more variable denudation rates through time, depending on the residence time of the landslide regolith and the frequency of earthquakes, than shallower basins. Landscapes with more variation in basin relief could enhance the temporal perturbations in denudation rates produced by earthquakes, but the contribution of tectonics to long-term variation is difficult to isolate.

We also explore the averaging time required to reach the mean exhumation rate in model runs with different contributions from seismic and aseismic uplift. After a large earthquake that produces tens of centimetres of instantaneous weathering of SB, exhumation rates measured with different averaging times converge to the long-term mean rate within hundreds to thousands of years (Fig. 5a). Bedrock surface exhumation rates are impacted by both surface uplift and lowering rates, so as a result the timescale of the perturbation is impacted by the dominant form of uplift in the mountain range. As a result, the more dominant coseismic uplift is in a landscape, the longer the recording time needs to be before a reliable exhumation record can be made. Landscapes with more frequent earthquakes have more variable exhumation rates which require longer averaging times to achieve accurate measurements of the long-term average exhumation rate. Regardless of the frequency of earthquakes in a mountain range, the events of the greatest magnitude remain uncommon while being the dominant contributors to weathering. Hence, the relationship between exhumation rate and averaging time is consistent with the Sadler effect that has been described for sedimentary systems (Schumer and Jerolmack, 2009). Unlike sedimentary systems, where it is possible to measure sedimentation rates from timescales of seconds to millions of years, there are few measures of exhumation that we can use and many of these have long averaging times. Thermochronometric methods average across timescales of 105–107 years, much longer than the recurrence times of individual earthquakes. There is a possibility that cosmogenic radionuclide analysis records individual earthquakes, where earthquake-triggered landslides weather bedrock that has a low cosmogenic concentration (Wang et al., 2017; West et al., 2014), although enhanced erosion during and immediately after an earthquake complicates the cosmogenic signal in practice.

Despite representing close to half of the weathering flux of mountain belts, stochastic earthquakes still remain substantively missing from our models of the development of mountain belts. The modelling here demonstrates that stochastic uplift and exhumation by large earthquakes is likely to be averaged away across the timescale of most thermochronometers, with the variability in uplift and exhumation representing around 15 % of the average exhumation rate. Even when a large earthquake occurs at the edge of a mountain belt, as has occurred in the Wenchuan region, the variable exhumation signal is further shredded by sediment transport by floods and debris flows, such that even sinks that are within 40 km of the epicentre show limited evidence for large earthquakes (Zhang et al., 2019). This result along with our model demonstrates the importance of understanding the processes by which landslide sediment are mobilised out of catchments and the time taken for these processes. Without improved understanding of the cascading nature of sediment transport from catchments, it is unlikely we will be able to identify earthquakes other than within smaller basins or sinks immediately adjacent to the epicentre (Howarth et al., 2012). Large basins (greater than 10 000 km2) have been shown to be large enough to average out the perturbations in cosmogenic radionuclide concentrations caused by large landsliding events (Dingle et al., 2018; Marc et al., 2019). While the largest basins are able to offer reliable estimates of long-term erosion rates, the detrital cosmogenic nuclide concentrations from smaller basins will be more affected by bedrock landsliding caused by earthquakes. Therefore, smaller basins could be suitable targets to recognise variations in erosion rates due to earthquakes. Our model suggests that the impact of a large earthquake is not necessarily big enough to perturb the denudation rate of an orogen for the whole of a cosmogenic nuclide recording time. The combination of averaging times, shredding, and the relative contributions of large earthquakes to long-term exhumation rates helps us to understand the lack of clear signatures for single earthquakes in sedimentary or exhumation records.

4 Conclusions and implications

Our simulations show that the regolith generated by large earthquakes can reduce the rate of weathering and exhumation of rock due to its potentially long residence time on hillslopes. Reducing exhumation rates also increases the uplift rate of the bedrock surface, but these effects are small when compared to the role of the magnitude and frequency of earthquakes. These results demonstrate that background tectonic processes and rates are the dominant control on surface uplift, while the more important role for large earthquakes is their control on weathering. Small earthquakes contribute very little to both uplift and weathering resulting in below-average rock exhumation being recorded if a large earthquake does not occur during the averaging time of the exhumation record. While large earthquakes produce higher-than-average rock exhumation rates, the slow removal of regolith from the orogen reduces the magnitude and timescale of the signal. The relatively long timescales of exhumation records prevent the recording of orogen-scale variation in exhumation due to a single earthquake. A better understanding of the controls on bedrock weathering by earthquake-triggered landslides is required to identify signals of earthquakes in sedimentary records. Higher-resolution exhumation records and the growing recognition of the complex nature of exporting landslide sediment from mountainous catchments will help to explore this problem.

Data availability

All data for Figs. 2–5 are generated via the model described in this paper, the original code of which is available upon request. The data from Fig. 6 can be found in the Harel et al. (2016) paper.

Author contributions

All authors contributed to the writing and ideas of the paper. The model and project design were produced with discussion from all authors. ORF wrote the model with assistance from DEJH, TCH, and AJH. TCH, XF, and RH formulated the overall project aims.

Competing interests

The authors declare that they have no conflict of interest.


Daniel E. J. Hobley was supported in part by an EU Marie Skłodowska-Curie Actions Sêr Cymru COFUND fellowship. The manuscript has benefited greatly from thorough comments of Alex Densmore and an anonymous reviewer.

Financial support

This research has been supported by the Newton Fund (grant no. NE/N012240/1), the Natural Environmental Research Council (grant no. NE/N012240/1), and the Economic and Social Research Council (grant no. NE/N012240/1).

Review statement

This paper was edited by Susan Conway and reviewed by Alexander Densmore and one anonymous referee.


Andermann, C., Crave, A., Gloaguen, R., Davy, P., and Bonnet, S.: Connecting source and transport: Suspended sediments in the Nepal Himalayas, Earth Planet. Sc. Lett., 351–352, 158–170,, 2012. 

Avouac, J. P.: Dynamic Processes in Extensional and Compressional Settings – Mountain Building: From Earthquakes to Geological Deformation, Treat. Geophys., 6, 377–439, 2007. 

Balco, G., Stone, J. O., Lifton, N. A., and Dunai, T. J.: A complete and easily accessible means of calculating surface exposure ages or erosion rates from 10Be and 26Al measurements, Quatern. Geochronol., 3, 174–195,, 2008. 

Bennett, G. L., Molnar, P., McArdell, B. W., and Burlando, P.: A probabilistic sediment cascade model of sediment transfer in the Illgraben, Water Resour. Res., 50, 1225–1244,, 2014. 

Blöthe, J. H. and Korup, O.: Millennial lag times in the Himalayan sediment routing system, Earth Planet. Sc. Lett., 382, 38–46,, 2013. 

Bonilla, M. G., Mark, R. K., and Lienkaemper, J. J.: Statistical relations among earthquake magnitude, surface rupture length and surface fault displacement, B. Seismol. Soc. Am., 74, 2379–2411, 1984. 

Braucher, R., Merchel, S., Borgomano, J., and Bourlès, D. L.: Production of cosmogenic radionuclides at great depth: A multi element approach, Earth Planet. Sc. Lett., 309, 1–9,, 2011. 

Clark, M. K., Bush, J. W. M., and Royden, L. H.: Dynamic topography produced by lower crustal flow against rheological strength heterogeneities bordering the Tibetan Plateau, Geophys. J. Int., 162, 575–590,, 2005. 

Croissant, T., Lague, D., Steer, P., and Davy, P.: Rapid post-seismic landslide evacuation boosted by dynamic river width, Nat. Geosci., 1, 680–684,, 2017. 

Croissant, T., Steer, P., Lague, D., Davy, P., Jeandet, L., and Hilton, R. G.: Seismic cycles, earthquakes, landslides and sediment fluxes: Linking tectonics to surface processes using a reduced-complexity model, Geomorphology, 339, 87–103,, 2019. 

Dai, F. C., Xu, C., Yao, X., Xu, L., Tu, X. B., and Gong, Q. M.: Spatial distribution of landslides triggered by the 2008 Ms 8.0 Wenchuan earthquake, China, J. Asian Earth Sci., 40, 883–895,, 2011. 

Densmore, A. L., Ellis, M. A., Li, Y., Zhou, R., Hancock, G. S., and Richardson, N.: Active tectonics of the Beichuan and Pengguan faults at the eastern margin of the Tibetan Plateau, Tectonics, 26, 1–17,, 2007. 

Densmore, A. L., Li, Y., Richardson, N. J., Zhou, R., Ellis, M., and Zhang, Y.: The role of late quaternary upper-crustal faults in the 12 may 2008 Wenchuan earthquake, B. Seismol. Soc. Am., 100, 2700–2712,, 2010. 

Densmore, A. L., Parker, R. N., Rosser, N. J., De Michele, M., Yong, L., Runqiu, H., Whadcoat, S., and Petley, D. N.: Reply to “Isostasy can't be ignored”, Nat. Geosci., 5, 83–84,, 2012. 

Dingle, E. H., Sinclair, H. D., Attal, M., Rodés, Á., and Singh, V.: Temporal variability in detrital 10Be concentrations in a large Himalayan catchment, Earth Surf. Dynam., 6, 611–635,, 2018. 

England, P. and Molnar, P.: Surface uplift, uplift of rocks, and exhumation of rocks, Geology, 18, 1173–1177,<1173:SUUORA>2.3.CO, 1990. 

Fan, X., van Westen, C. J., Korup, O., Gorum, T., Xu, Q., Dai, F., Huang, R., and Wang, G.: Transient water and sediment storage of the decaying landslide dams induced by the 2008 Wenchuan earthquake, China, Geomorphology, 171–172, 58–68,, 2012. 

Fan, X., Domènech, G., Scaringi, G., Huang, R., Xu, Q., Hales, T. C., Dai, L., Yang, Q., and Francis, O.: Spatio-temporal evolution of mass wasting after the 2008 Mw7.9 Wenchuan Earthquake revealed by a detailed multi-temporal inventory, Landslides, 15, 2325–2341,, 2018. 

Giardini, D., Grunthal, G., Shedlock, K. M., and Peizhen, Z.: The GSHAP Global Seismic Hazard Map, Ann. Di Geofis., 42, 1225–1230, 1999. 

Gosse, J. and Philips, F.: Terrestrial in situ cosmogenic nuclides:theory and application, Quaternary Sci. Rev., 20, 1475–1560, 2001. 

Granger, D. E. and Muzikar, P. F.: Dating sediment burial with in situ-produced cosmogenic nuclides: Theory, techniques, and limitations, Earth Planet. Sc. Lett., 188, 269–281,, 2001. 

Hales, T. C., Abt, D. L., Humphreys, E. D., and Roering, J. J.: A lithospheric instability origin for Columbia River flood basalts and Wallowa Mountains uplift in northeast Oregon, Nature, 438, 842–845,, 2005. 

Harel, M., Mudd, S. M., and Attal, M.: Geomorphology Global analysis of the stream power law parameters based on worldwide Be denudation rates, Geomorphology, 268, 184–196,, 2016. 

Heimsath, A. M., Dietrichs, W. E., Nishiizuml, K., and Finkel, R. C.: The soil production function and landscape equilibrium, Nature, 388, 358–361,, 1997. 

Horton, A. J., Hales, T. C., Ouyang, C., and Fan, X.: Identifying post-earthquake debris flow hazard using Massflow, Eng. Geol., 258, 105134,, 2019. 

Hovius, N., Meunier, P., Lin, C. W., Chen, H., Chen, Y. G., Dadson, S., Horng, M. J., and Lines, M.: Prolonged seismically induced erosion and the mass balance of a large earthquake, Earth Planet. Sc. Lett., 304, 347–355,, 2011. 

Howarth, J. D., Fitzsimons, S. J., Norris, R. J., and Jacobsen, G. E.: Lake sediments record cycles of sediment flux driven by large earthquakes on the Alpine fault, New Zealand, Geology, 40, 1091–1094,, 2012. 

Hubbard, J. and Shaw, J. H.: Uplift of the Longmen Shan and Tibetan plateau, and the 2008 Wenchuan (M=7.9) earthquake, Nature, 458, 194–197,, 2009. 

Iverson, R. M.: Landslide triggering by rain infiltration, Water Resour. Res., 36, 1897–1910,, 2000. 

Jerolmack, D. J. and Paola, C.: Shredding of environmental signals by sediment transport, Geophys. Res. Lett., 37, L19401,, 2010. 

Keefer, D. K.: The importance of earthquake-induced landslides to long-term slope erosion and slope-failure hazards in seismically active regions, Geomorphology, 10, 265–284,, 1994. 

Kirby, E., Whipple, K. X., Tang, W., and Chen, Z.: Distribution of active rock uplift along the eastern margin of the Tibetan Plateau: Inferences from bedrock channel longitudinal profiles, J. Geophys. Res.-Sol. Ea., 108, 2217,, 2003. 

Kirchner, J. W., Finkel, R. C., Riebe, C. S., Granger, D. E., Clayton, J. L., King, J. G., and Megahan, W. F.: Mountain erosion over 10 yr, 10 k.y., and 10 m.y. time scales, Geology, 29, 591–594,<0591:MEOYKY>2.0.CO;2, 2001. 

Korup, O. and Clague, J. J.: Natural hazards, extreme events, and mountain topography, Quaternary Sci. Rev., 28, 977–990,, 2009. 

Larsen, I. J., Montgomery, D. R., and Korup, O.: Landslide erosion controlled by hillslope material, Nat. Geosci., 3, 247–251,, 2010. 

Li, G., West, A. J., Densmore, A. L., Jin, Z., Parker, R. N., and Hilton, R. G.: Seismic mountain building: Landslides associated with the 2008 Wenchuan earthquake in the context of a generalized model for earthquake volume balance, Geochem. Geophy. Geosy., 15, 833–844,, 2014. 

Li, G., West, A. J., Densmore, A. L., Hammond, D. E., Jin, Z., Zhang, F., Wang, J., and Hilton, R. G.: Connectivity of earthquake-triggered landslides with the fluvial network: Implications for landslide sediment transport after the 2008 Wenchuan earthquake, J. Geophys. Res.-Earth, 121, 703–724,, 2016. 

Li, G., West, A. J., Densmore, A. L., Jin, Z., Zhang, F., Wang, J., Clark, M., and Hilton, R. G.: Earthquakes drive focused denudation along a tectonically active mountain front, Earth Planet. Sc. Lett., 472, 253–265,, 2017. 

Li, G., West, A. J., and Qiu, H.: Competing effects of mountain uplift and landslide erosion over earthquake cycles, J. Geophys. Res.-Sol. Ea., 124, 5101–5133,, 2019. 

Li, Z., Liu-Zeng, J., Almeida, R., Hubbard, J., Sun, C., and Yi, G.: Re-evaluating seismic hazard along the southern Longmen Shan, China: Insights from the 1970 Dayi and 2013 Lushan earthquakes, Tectonophysics, 717, 519–530,, 2017. 

Liu-Zeng, J., Zhang, Z., Wen, L., Tapponnier, P., Sun, J., Xing, X., Hu, G., Xu, Q., Zeng, L., Ding, L., Ji, C., Hudnut, K. W., and van der Woerd, J.: Co-seismic ruptures of the 12 May 2008, Ms8.0 Wenchuan earthquake, Sichuan: East-west crustal shortening on oblique, parallel thrusts along the eastern edge of Tibet, Earth Planet. Sc. Lett., 286, 355–370,, 2009. 

Malamud, B. D., Turcotte, D. L., Guzzetti, F., and Reichenbach, P.: Landslide inventories and their statistical properties, Earth Surf. Process. Land., 29, 687–711,, 2004a. 

Malamud, B. D., Turcotte, D. L., Guzzetti, F., and Reichenbach, P.: Landslides, earthquakes, and erosion, Earth Planet. Sc. Lett., 229, 45–59,, 2004b. 

Marc, O., Hovius, N., Meunier, P., Uchida, T., and Hayashi, S.: Transient changes of landslide rates after earthquakes, Geology, 43, 883–886,, 2015. 

Marc, O., Hovius, N., Meunier, P., Gorum, T., and Uchida, T.: A seismologically consistent expression for the total area and volume of earthquake-triggered landsliding, J. Geophys. Res.-Earth, 121, 640–663,, 2016a. 

Marc, O., Hovius, N., and Meunier, P.: The mass balance of earthquakes and earthquake sequences, Geophys. Res. Lett., 43, 3708–3716,, 2016b. 

Marc, O., Behling, R., Andermann, C., Turowski, J. M., Illien, L., Roessner, S., and Hovius, N.: Long-term erosion of the Nepal Himalayas by bedrock landsliding: the role of monsoons, earthquakes and giant landslides, Earth Surf. Dynam., 7, 107–128,, 2019. 

Meade, B. J.: The signature of an unbalanced earthquake cycle in Himalayan topography?, Geology, 38, 987–990,, 2010. 

Meunier, P., Uchida, T., and Hovius, N.: Landslide patterns reveal the sources of large earthquakes, Earth Planet. Sc. Lett., 363, 27–33,, 2013. 

Molnar, P.: Isostasy can't be ignored, Nat. Geosci., 5, p. 83,, 2012. 

Molnar, P., England, P., and Martinod, J.: Mantle dynamics, uplift of the Tibetan plateau, and the Indean monsoon, Rev. Geophys., 31, 357–396, 1993. 

Molnar, P., England, P. C., and Jones, C. H.: Mantle dynamics, isostasy, and the support of high terrain, J. Geophys. Res.-Sol. Ea., 120, 1932–1957,, 2015. 

Montgomery, D. R. and Brandon, M. T.: Topographic controls on erosion rates in tectonically active mountain ranges, Earth Planet. Sc. Lett., 201, 481–489,, 2002. 

Niemi, N. A., Oskin, M., Burbank, D. W., Heimsath, A. M., and Gabet, E. J.: Effects of bedrock landslides on cosmogenically determined erosion rates, Earth Planet. Sc. Lett., 237, 480–498,, 2005. 

Ouimet, W. B.: Landslides associated with the May 12, 2008 Wenchuan earthquake: Implications for the erosion and tectonic evolution of the Longmen Shan, Tectonophysics, 491, 244–252,, 2010. 

Ouimet, W. B., Whipple, K. X., Royden, L. H., Sun, Z., and Chen, Z.: The influence of large landslides on river incision in a transient landscape: Eastern margin of the Tibetan Plateau (Sichuan, China), B. Geol. Soc. Am., 119, 1462–1476,, 2007. 

Ouimet, W. B., Whipple, K. X., and Granger, D. E.: Beyond threshold hillslopes: Channel adjustment to base-level fall in tectonically active mountain ranges, Geology, 37, 579–582,, 2009. 

Parker, R. N., Densmore, A. L., Rosser, N. J., De Michele, M., Li, Y., Huang, R., Whadcoat, S., and Petley, D. N.: Mass wasting triggered by the 2008 Wenchuan earthquake is greater than orogenic growth, Nat. Geosci., 4, 449–452,, 2011. 

Parker, R. N., Hancox, G. T., Petley, D. N., Massey, C. I., Densmore, A. L., and Rosser, N. J.: Spatial distributions of earthquake-induced landslides and hillslope preconditioning in the northwest South Island, New Zealand, Earth Surf. Dynam., 3, 501–525,, 2015. 

Pearce, A. J. and Watson, A. J.: Effects of earthquake-induced landslides on sediment budget and transport over a 50-yr period, Geology, 14, 52–55, 1986. 

Robinson, T. R., Davies, T. R. H., Wilson, T. M., and Orchiston, C.: Coseismic landsliding estimates for an Alpine Fault earthquake and the consequences for erosion of the Southern Alps, New Zealand, Geomorphology, 263, 71–86,, 2016. 

Royden, L. H., Burchfiel, B. C., King, R. W., Chen, Z., Shen, F., and Liu, Y.: Surface deformation and lower crust flow in eastern Tibet, Science, 276, 788–791, 1997. 

Schumer, R. and Jerolmack, D. J.: Real and apparent changes in sediment deposition rates through time, J. Geophys. Res.-Sol. Ea., 114, 1–12,, 2009. 

Simpson, G.: Accumulation of permanent deformation during earthquake cycles on reverse faults, J. Geophys. Res.-Sol. Ea., 120, 1958–1974,, 2015. 

Stark, C. P., Hovius, N., and Stark, C. P.: The characterization of landslide size distributions, Geophys. Res. Lett., 28, 1091–1094, 2001. 

Stolle, A., Bernhardt, A., Schwanghart, W., Hoelzmann, P., Adhikari, B. R., Fort, M., and Korup, O.: Catastrophic valley fills record large Himalayan earthquakes, Pokhara, Nepal, Quaternary Sci. Rev., 177, 88–103,, 2017. 

Turcotte, D. L. and Schubert, G.: Geodynamics, 2nd Edn., Cambridge University Press, Cambridge, 2002.  

Valagussa, A., Marc, O., Frattini, P., and Crosta, G. B.: Seismic and geological controls on earthquake-induced landslide size, Earth Planet. Sc. Lett., 506, 268–281,, 2019. 

Wang, W., Godard, V., Liu-Zeng, J., Scherler, D., Xu, C., Zhang, J., Xie, K., Bellier, O., Ansberque, C., and de Sigoyer, J.: Perturbation of fluvial sediment fluxes following the 2008 Wenchuan earthquake, Earth Surf. Proc. Land., 42, 2611–2622,, 2017. 

Wells, D. L. and Coppersmith, K. J.: New empical relationship between magnitude, rupture length, rupture width, rupture area, and surface displacement, B. Seismol. Soc. Am., 84, 974–1002, 1994. 

West, A. J., Hetzel, R., Li, G., Jin, Z., Zhang, F., Hilton, R. G., and Densmore, A. L.: Dilution of 10Be in detrital quartz by earthquake-induced landslides: Implications for determining denudation rates and potential to provide insights into landslide sediment dynamics, Earth Planet. Sc. Lett., 396, 143–153,, 2014. 

Yanites, B. J., Tucker, G. E., and Anderson, R. S.: Numerical and analytical models of cosmogenic radionuclide dynamics in landslide-dominated drainage basins, J. Geophys. Res.-Earth, 114, F01007,, 2009. 

Yanites, B. J., Tucker, G. E., Mueller, K. J., and Chen, Y.-G.: How rivers react to large earthquakes: Evidence from central Taiwan, Geology, 38, 639–642,, 2010. 

Zhang, F., Jin, Z., West, A. J., An, Z., Hilton, R. G., Wang, J., Li, G., Densmore, A. L., Yu, J., Qiang, X., Sun, Y., Li, L., Gou, L., Xu, Y., Xu, X., Liu, X., Pan, Y., and You, C.-F.: Monsoonal control on a delayed response of sedimentation to the 2008 Wenchuan earthquake, Sci. Adv., 5, eaav7110,, 2019. 

Zhang, S., Zhang, L., Lacasse, S., and Nadim, F.: Evolution of Mass Movements near Epicentre of Wenchuan Earthquake, the First Eight Years, Sci. Rep., 6, 1–9,, 2016. 

Short summary
Large earthquakes can build mountains by uplifting bedrock, but they also erode them by triggering large volumes of coseismic landsliding. Using a zero-dimensional numerical model, we identify that the storage of sediment produced by earthquakes can affect surface uplift and exhumation rates across the mountain range. However, the storage also reduces the time span at which the impact of the earthquake can be measured, preventing the recognition of single earthquakes in many long-term records.