Quantifying retrogressive thaw slump mass wasting and carbon mobilisation on the Qinghai-Tibet Plateau using multi-modal remote sensing

Retrogressive Thaw Slumps (RTS) are slope failures triggered by permafrost thaw that occur in ground-ice-rich regions of the Arctic and the Qinghai-Tibet Plateau (QTP). A strong warming trend has amplified RTS activity on the QTP in recent years. Although the region currently acts as a carbon sink, its permafrost-covered area (40 %) contains substantial soil organic carbon (SOC) stocks. Intensifying thaw-driven mass wasting may transform the QTP into a net carbon source by mobilising previously frozen SOC and enhancing decomposition. Yet, regional remote sensing studies have not yet quantified RTS mass wasting, including material erosion volumes and associated SOC mobilisation. Analysing time-series data from digital elevation models (DEMs) enables direct observation of RTS activity by measuring changes in active area, eroded material volume, and the overall magnitude of surface change. However, most available DEM sources lack the spatial resolution and temporal frequency required for comprehensive RTS monitoring. In contrast, optical data provide higher spatial resolution and more frequent observations, but lack elevation information. Here, we evaluated RTS mass wasting across the QTP from 2011 to 2020 by combining DEMs derived from bistatic Interferometric Synthetic Aperture Radar (InSAR) observations of the TanDEM-X mission with annual RTS inventories generated from high-resolution optical satellite imagery and geophysical soil property data to estimate erosion volume, ground ice loss, and SOC mobilisation. We estimated that RTS activity on the QTP during 2011–2020 relocated 5.02${}_{\mathrm{0.75}}^{\mathrm{25.35}}\times {\mathrm{10}}^{\mathrm{7}}$ m³ previously frozen material, resulting in a loss of 3.58${}_{\mathrm{0.28}}^{\mathrm{28.20}}\times {\mathrm{10}}^{\mathrm{6}}$ m³ of ground ice, and mobilised 2.78${}_{\mathrm{0.11}}^{\mathrm{7.98}}\times {\mathrm{10}}^{\mathrm{8}}$ kg C of organic carbon. We found a reliable power-law scaling relationship between RTS area in the optical RTS inventory and calculated volume change, with α values ranging from 1.20±0.01 to 1.30±0.01 (R²=0.87, p<0.001) depending on the regression model used, which may readily transform planimetric RTS area into volume estimates at scale on the QTP. Despite the relatively recent initiation and smaller size of RTSs on the QTP, material erosion and SOC mobilisation over the past decade exceeded levels in some Siberian Arctic regions, but remained up to 10 times lower than hotspots in the Canadian High Arctic. While current RTS impacts on the QTP are relatively modest, affecting < 0.01 % of the total permafrost area and contributing approximately 0.1 % to the regional carbon budget, the accelerating rates of RTS activity indicate that this phenomenon could become increasingly significant in the future. Our findings highlight the importance of regional studies in advancing our understanding of permafrost thaw-driven changes to the carbon dynamics of rapidly changing permafrost ecosystems.

1 Introduction

Permafrost regions are rapidly warming, causing widespread thaw and degradation (Biskaborn et al., 2019; Farquharson et al., 2019; O'Neill et al., 2023). The permafrost thaw is triggered by a long-term increase in air temperatures (Smith et al., 2022) and is further amplified by various short-term disturbances, such as periods of extreme summer temperatures, high-intensity rainfall, hydrological changes, wildfires, and anthropogenic impacts (Grosse et al., 2011; Hjort et al., 2022; Holloway et al., 2020; Bernhard et al., 2022b; Kokelj et al., 2015). Although approximately half of the global soil organic carbon (SOC) stock is stored in permafrost soils in the Northern Hemisphere (Mishra et al., 2021; Schuur et al., 2022), ongoing permafrost decline is expected to accelerate SOC decomposition and greenhouse gas emissions, potentially triggering significant climate feedbacks (Schuur et al., 2022; Yi et al., 2025). Due to the lack of large-scale observations and the complexity of permafrost thaw processes, climate models do not account for the potential of such processes when simulating permafrost carbon feedbacks (Yi et al., 2025; Schuur et al., 2015, 2022). Existing Earth system models exhibit significant limitations in the accounting of soil organic carbon (SOC) and in the prediction of future changes for global permafrost regions (Turetsky et al., 2020; Virkkala et al., 2021; Nitzbon et al., 2020).

The largest high-altitude permafrost zone is the Qinghai-Tibet Plateau (QTP) with a total extent of 1.06×10⁶ km² at mean elevations greater than 4000 m (Wang and French, 1994; Liu and Chen, 2000; Zou et al., 2017). Similarly to high-latitude permafrost regions, the QTP is one of the most climate-sensitive regions on Earth (Liu and Chen, 2000) and has experienced a pronounced warming trend in recent decades, with an average increase in air temperatures of 0.035 °C a⁻¹ (Yao et al., 2019). The warming trend affects the thermal state of the permafrost: the active layer thickness (ALT) and ground temperatures have increased, while the regional permafrost extent declined (Cheng and Wu, 2007; Wu and Zhang, 2008; Zhao et al., 2021; Ran et al., 2022). Hence, the QTP is susceptible to permafrost thaw processes that have a substantial impact on the environment and communities, including threats to local transport and energy infrastructure, ecosystems and hydrology, as well as regional carbon budgets and water storage capacities (Luo et al., 2019; Li et al., 2022; Mu et al., 2017; Zhao et al., 2020; Yi et al., 2025; Chen et al., 2024b). Permafrost thaw together with the large SOC stocks make the QTP a potentially considerable carbon source and an important region to monitor permafrost thaw processes (Ran et al., 2022; Chen et al., 2024a; Yi et al., 2025).

Retrogressive Thaw Slumps (RTS) are permafrost landforms that occur in ice-rich permafrost terrain when ground ice is exposed, allowing for rapid thaw and downslope movement of the resulting debris (Burn and Lewkowicz, 1990; Kokelj and Jorgenson, 2013; Nesterova et al., 2024; Harris, 1988; CPA et al., 2024). The landform can expand successively upslope with time due to continuous exposure and thawing of massive and segregated ground ice, thus eroding steep headwalls and mobilising thawed material downslope and potentially into nearby streams and rivers (Nesterova et al., 2024; Kokelj et al., 2021). A gradual decline in headwall height with upslope growth and material accumulation can also result in stabilisation and recovery (Burn and Lewkowicz, 1990; van der Sluijs et al., 2023). RTS can reactivate in complex polycyclic ways (Krautblatter et al., 2024). Climate warming and human disturbance have intensified RTS activity not only in the Arctic (Lantz and Kokelj, 2008; Bernhard et al., 2022a; van der Sluijs et al., 2023), but also on the QTP. In recent years, the plateau experienced strong expansion and initiation rates of RTS mainly on gentle north-facing slopes with fine-grained soils and high ground ice content. More than 30 % of all RTS activity is observed in the Beiluhe River Basin located in the central QTP, where most activity started after 2010 (Luo et al., 2019, 2022; Huang et al., 2020; Xia et al., 2022, 2024). RTS retreat rates are relatively high, with mean rates up to 25 m a⁻¹ (2017–2019) (Huang et al., 2021), although these are similar to other highly active RTS sites in Alaska, northwest Canada, the Canadian High Arctic, and Siberia (Ward Jones et al., 2019; Lacelle et al., 2015, 2010; Lantuit and Pollard, 2008; Lantz and Kokelj, 2008; Leibman et al., 2021).

Due to their complex spatiotemporal dynamics, monitoring RTS activity and assessing their impact on regional carbon cycling remains challenging and is still associated with considerable uncertainties. Most studies have focused on thermokarst hotspots in the central and northeast regions of the QTP, conducting local to subregional analyses using high-resolution satellite imagery and (semi)automatic detection algorithms (Luo et al., 2022; Huang et al., 2021; Xia et al., 2022). A recent study provided an annual inventory of more than 3000 RTS features across the QTP from 2016 to 2022, using high-resolution optical PlanetScope imagery and a semiautomatic detection approach. This effort represents the first high-quality regional-scale dataset that describes RTS initiation and planimetric expansion (Xia et al., 2024). Nevertheless, to quantify RTS-induced mass wasting and evaluate the potential implications on permafrost carbon mobilisation, additional datasets are required - particularly those capturing lateral and vertical change and soil properties. By deriving elevation change from pairs of high-resolution digital elevation models (DEMs) with a temporal baseline within the affected planimetric RTS boundary, volume changes of the eroded material can be obtained by air- or spaceborne LiDAR, stereo optical photogrammetry, or bistatic Interferometric Synthetic Aperture Radar (InSAR) measurements (Lantuit and Pollard, 2005; Van Der Sluijs et al., 2018; Ramage et al., 2018; Dai et al., 2024; Bernhard et al., 2020). Data from the German TanDEM-X mission's X-band bistatic SAR enables the creation of global multitemporal DEMs with approximately 10 m spatial resolution and 2–3 m vertical accuracy (Krieger et al., 2007) and has been shown to be suitable for pan-Arctic monitoring of RTS mass wasting (Bernhard et al., 2022b; Maier et al., 2025). Using in-situ measurements or modelled estimates of SOC stocks and ground ice content, we can combine this information with volumetric change estimates to estimate the amount of formerly frozen organic carbon that has been mobilised. Such approaches have been applied to assess carbon mobilisation due to coastal erosion and slumping in Canada (Ramage et al., 2018), a severe heat wave in Siberia (Bernhard et al., 2022b), and more recently at large thaw-driven mass wasting sites across the pan-Arctic (Dai et al., 2025). However, to date, there have been no regional-scale empirical estimates of either material erosion volume or SOC mobilisation for the QTP.

The volume of eroded material of an RTS scales with its area following power-law relations that characterise its growth dynamics. Several studies adapted the so-called area-volume or allometric scaling from temperate landslide research (Jaboyedoff et al., 2020). Commonly, area-volume scaling is performed using an ordinary least squares (OLS) approach to fit a linear model to the log-transformed RTS area and volume (change) to obtain scaling coefficients. However, distinct differences in the estimated scaling laws can be present between geographic regions and based on the methodological approach used. Bernhard et al. (2022a) used an orthogonal distance regression (ODR) (Boggs and Rogers, 1989) instead of OLS to fit the straight line to the log-transformed RTS area and volume change based on TanDEM-X DEM pairs (2010–2016) assuming that both volume change δV and area change δA are affected by measurement error. The authors report a scaling coefficient of 1.15 for several Arctic sites. Dai et al. (2025) reported a pan-Arctic scaling coefficient of 1.30 based on OLS and ArcticDEM pairs (2012–2022) while Kokelj et al. (2021) and van der Sluijs et al. (2023) report coefficients of 1.36 and 1.41, respectively, in the Canadian Arctic based on OLS and pre-disturbance terrain reconstruction (until 2018). A recent study estimated a scaling coefficient of 1.20 for almost 1500 RTS on the QTP (Ma et al., 2025) based on DEM mosaics and commercial stereo-optical DEMs with varying dates (until 2021–2025). Robust empirical scaling relationships can be helpful to potentially constrain regional-scale estimates on material erosion and carbon mobilisation when only optically derived RTS area estimates are available. However, differences in the scaling model and temporal and spatial resolution of the elevation data impair the (inter-)regional transferability of the estimated coefficients.

In this study, we present the first regional empirical analysis on RTS mass wasting due to RTS activity on the QTP during the last decade. The elevation change from TanDEM-X-derived DEMs between 2011 and 2020, in combination with the high-quality RTS delineations of Xia et al. (2024) and soil property datasets for the QTP, allowed us to estimate the volume of eroded material and the related allometric scaling coefficients, associated ground ice loss, and mobilisation of SOC. We aim to show that the combination of multimodal and multitemporal datasets allows for a more detailed analysis of RTS mass wasting dynamics and further increases our understanding of the regional carbon budget impacts.

2 Data and methods

2.1 Study site

Our study region, the QTP, is located between 26 and 38° N in the south-west of China at average elevations higher than 4000 m above sea level (Fig. 1a). Permafrost covers 40 % of the plateau (Wang and French, 1994; Liu and Chen, 2000; Zou et al., 2017). The permafrost ground ice content averages around 30 %, decreasing spatially from north to south and west to east (Zou et al., 2024) (Fig. 1b). Compared to the Arctic, the ALT is high ($\stackrel{\mathrm{‾}}{\mathrm{ALT}}=\mathrm{2.34}$ m) (Ran et al., 2022), while permafrost thickness is relatively low (<60–350 m) (Zhao et al., 2020). A dry and cold climate in the northwest transitions to a warmer and wetter climate in the southeast of the plateau (Chen et al., 2015). The QTP permafrost also stores large amounts of SOC with a median estimate of 1.41×10¹³ kg C (or 14.1 Pg C) for the top 3 m and 4.92×10¹³ kg C (or 49.2 Pg C) for the upper 25 m of soils (Wang et al., 2020; Chen et al., 2024a). The SOC stocks increase from west to east and from north to south (Wang et al., 2021; Chen et al., 2024a) (Fig. 1c).

Figure 1 Overview of the study region, the Qinghai-Tibet Plateau (QTP). All data is gridded into 50×50 km² units for visualisation purposes. (a) Distribution and spatial density of RTS (Xia et al., 2024). We divided the QTP into five subregions (West, West-Central, Central, East, Northeast) and distributed validation sites (400 km²) across the QTP. The background depicts the spatial distribution of permafrost terrain on the QTP (Zou et al., 2017). (b) The ground ice content on the QTP for an exemplary depth layer (2–3 m) is highest in the central part of the plateau (Zou et al., 2024). (c) Soil organic carbon (SOC) stocks aggregated for the first 3 m depth show an increasing trend from the northwest to the southeast QTP (Wang et al., 2021).

Based on the spatial clustering of RTS identified by Xia et al. (2024), we divided the study area into five subregions, including: West, West-Central, Central, East, and Northeast (Fig. 1a). These subregions are used to analyse the spatial patterns of material erosion and mobilisation of SOC induced by RTS activity, and area-volume-scaling across the QTP. To validate the estimated volumes of RTS material erosion derived from integrated optical and elevation remote sensing data, we established the following validation sites covering 400 km² each (Fig. 1a): Western Kunlun, Gaize, Southern Nima, Beiluhe River Basin, and Qilian Mountain. We selected sites based on (1) the presence of RTS activity located within diverse geographic and terrain conditions identified by Xia et al. (2024), (2) the availability of field observations, and (3) the sufficient coverage of TanDEM-X observations to ensure consistent data quality for validation purposes.

2.2 Workflow and data processing

We used bistatic TanDEM-X SAR observations to generate multi-temporal DEMs that cover all RTS locations in the RTS inventory of Xia et al. (2024). By differencing DEMs from 2011 and 2020 and combining the resulting elevation change with the annual high-resolution RTS inventory, we estimated the volume change of the eroded material induced by RTS activity on the QTP (Fig. 2a, b). Based on recent datasets on the permafrost state and soil conditions including ALT, ground ice content (GI) content, and SOC stocks, we modelled the material erosion volumes into annual SOC mobilisation rates for all RTS present on the QTP until 2020 (Fig. 2c). We evaluated the uncertainty in the estimated erosion volume and the derived properties, examining how spatial resolution affects errors in material erosion estimates (Fig. 2d). Similar to temperate landslides, scaling laws between the planimetric area and the erosion volume have been used to improve our understanding of the variability in geomorphology, process dynamics, and the drivers and controls of RTS.

Figure 2 Data processing and validation: (a) DEM generation from TanDEM-X SAR observations based on bistatic interferometric SAR (InSAR) processing and post-processing with data tiling, averaging, and differencing including coregistration routines resulting in elevation change maps between time T1 and T2. (b) We extract the elevation change (T2–T1) within the boundaries of the temporally matching optical RTS labels of Xia et al. (2024) divided into negative elevation change (ablation) and positive elevation change (deposition of thawed material) within the RTS delineation. (c) Estimation of the properties of RTS mass wasting including material erosion calculated from negative elevation change within the RTS boundary, ground ice (GI) loss, and SOC mobilisation. (d) Validation at selected sites by (1) evaluating how accurately we can estimate the active erosion area δA and the material erosion volume δV based on the optical RTS delineations and (2) by comparing approximate headwall heights obtained from very-high-resolution (VHR) DEMs of six RTS in the Beiluhe River Basin from photogrammetric drone surveys to our results.

2.2.1 Digital Elevation Model generation and processing

We used satellite observations from the German Synthetic Aperture Radar (SAR) mission TanDEM-X to generate temporally resolved DEMs based on bistatic SAR interferometry (InSAR) (Krieger et al., 2007; Bojarski et al., 2021). We excluded observations with height-of-ambiguity (HoA) values below 15 and above 80 m to guarantee a vertical accuracy between 2 and 3 m in flat areas (Martone et al., 2012; Bernhard et al., 2020). We used TanDEM-X observations throughout the year to achieve full spatial coverage on the QTP, although previous studies in the Arctic have used only winter data to avoid errors caused by dense and wet tundra vegetation during the growing season or melting snow packs (Bernhard et al., 2020, 2022a; Maier et al., 2025). The potential errors introduced by vegetation characteristics and snow cover are likely to be negligible in this study due to the commonly low canopy heights of the local vegetation (alpine meadows, arid desert, bare ground) (Wang et al., 2016; Xia et al., 2024), and shallow average snow depths well below the height sensitivity of the TanDEM-X DEMs (Che et al., 2008; Yang et al., 2020).

Using the global 12 m spatial resolution TanDEM-X DEM as a reference, pairs of bistatic SAR observations were processed with the GAMMA Remote Sensing software (Werner et al., 2000) to generate a series of DEM products following a standard InSAR processing workflow (Fig. 2a). Key steps include the generation of a differential interferogram, phase unwrapping, phase-to-height conversion, update of the reference DEM with the computed height difference, and the orthorectification or geocoding of the resulting DEM to map coordinates. Further information about the DEM generation pipeline can be found in Bernhard et al. (2020) and Maier et al. (2025). SAR shadow and layover areas, as well as the regions in the SAR images that experienced low coherence (<0.3), were not considered in subsequent processing (Bernhard et al., 2020; Maier et al., 2025). We estimated the associated random elevation error σ_h of each pixel for all generated DEMs using the interferometric coherence, HoA, and the multilook window of 4×4 pixels (Krieger et al., 2007; Rosen et al., 2000; Rodriguez and Martin, 1992). We reprojected all DEM products to a common horizontal coordinate system, WGS 84 / NSIDC EASE-Grid 2.0 Global (EPSG:6933) with an ellipsoidal vertical reference, resampled to 10 m spatial resolution, and split all data into 100 km²-tiles with a small spatial overlap to avoid any edge effects. We corrected the DEMs for vertical offsets and tilts and co-registered the DEM pairs using the Python package xDEM (Hugonnet et al., 2021). Due to the non-uniform temporal coverage of TanDEM-X observations on the QTP, we decided to calculate the elevation change δh between two time periods. We averaged all available DEMs for the time period T1 that spans one summer period (January 2011–April 2012) and for T2 that spans two summer seasons (September 2017–January 2020). In case of the existence of several DEMs in the same location during the same time period, we calculated a weighted average per pixel based on interferometric coherence before computing the elevation difference T2–T1 for each tile with RTS activity. For stable terrain, the resulting elevation change products were normally distributed around zero with a standard deviation representing the achievable vertical accuracy of the DEM pair. δh∼0 m indicates stable terrain, while negative and positive values reflect material ablation and deposition of thawed materials, respectively.

2.2.2 Multimodal RTS mass wasting quantification

Xia et al. (2024) created a high-quality regional inventory of annual RTS delineations between 2016 and 2022 based on a semi-automated deep learning approach with optical high-resolution PlanetScope imagery (Xia et al., 2022) (Fig. 2b). Spectral information in optical images distinguishes undisturbed from disturbed terrain, using differences in bare-ground and vegetation cover. RTS delineations derived from optical imagery often incorporate a broader area than the active ablation zone, including zones of recent activity and depositional sections of the slump floor (i.e., features not directly involved in ongoing material loss), while excluding stable zones of past disturbance masked by lush vegetation growth. Therefore, we calculated the volume of eroded material δV for each RTS by summing only the negative elevation changes δh_<0 m multiplied by the area of the DEM pixel (Fig. 2c).

$$\begin{array}{}\text{(1)}& \mathit{\delta }V=\sum \mathit{\delta }{h}_{<\mathrm{0}\mathrm{m}}\cdot \mathrm{100}{\mathrm{m}}^{\mathrm{2}}\end{array}$$

The DEMs of T2 may consist of up to two summer seasons of RTS activity, depending on the availability and spatial distribution of the observations. If the data from T2 consisted of observations from only one year, the optical RTS delineation from the matching year was selected as the RTS delineation for the mass wasting calculation. If observations from several years contributed to the DEM of T2, we assigned the optical RTS delineation that matched the most recent DEM to ensure we captured the full planimetric extent of the RTS. Nevertheless, if multiple observations contributed to the T2 DEM, the multi-year averaged erosion volume may underestimate the actual erosion volume.

Since previous studies on allometric scaling use different methods to fit a straight line to the log-transformed RTS area and volume change, we apply two common models: We use (1) an orthogonal distance regression (ODR) model (Boggs and Rogers, 1989) used by Bernhard et al. (2022a) for several North American and Siberian RTS sites, and (2) an ordinary least squares (OLS) approach applied by, for example, Kokelj et al. (2021), van der Sluijs et al. (2023), and Dai et al. (2025), to predict the eroded volume δV based on the planimetric area δA with an exponential scaling coefficient α and a scaling factor c (Jaboyedoff et al., 2020) for the time interval T1–T2:

$$\begin{array}{}\text{(2)}& \mathit{\delta }V=c\cdot \mathit{\delta }{A}^{\mathit{\alpha }}.\end{array}$$

In particular, δA stands for the RTS area that undergoes a negative elevation change or ablation within the monitoring period of T1–T2. To test the influence of different definitions of the RTS area between delineations on optical images and the elevation change maps, we also performed area-volume scaling with the entire RTS delineation area A_Xia, which does not only include zones of ablation δA but often also bare soils disturbed by mud flows and the deposition of thawed material in the slump floor (Fig. 2b).

To estimate how much SOC is stored in the previously frozen permafrost soils mobilised by RTS activity during the study period, ground conditions must be known or modelled (Fig. 2c). Therefore, we integrated existing datasets for the QTP that define (1) ALT (Ran et al., 2022), (2) GI content between 2 and 10 m (Zou et al., 2024), and (3) SOC stocks between 0 and 3 m depth (Wang et al., 2021). All datasets are sampled in 1 km²-cells. Ran et al. (2022) estimated permafrost thermal state variables, including ALT for the pan-Arctic permafrost region and the QTP with 452 field measurements with statistical learning models that achieve a root mean square error (RMSE) <1 m (1). Zou et al. (2024) estimated the GI content up to 10 m depth with a random forest algorithm based on climate, terrain and soil variables and 664 borehole records with R²>0.8 for all depth layers (Fig. 1b, (2)). Wang et al. (2021) applied a set of machine learning algorithms together with environmental variables and soil profile data from 572 field measurements to predict the spatial distribution of SOC in the upper 3 m, including an uncertainty layer on the QTP. The model achieved R² values between 0.66 for the upper 30 cm and 0.54 for the first metre (Fig. 1c, (3)). Since the SOC stocks dataset is limited to the upper 3 m, we applied an exponential decay model to extrapolate the values to lower depths where we propagated the uncertainty estimates of 2–3 m (Bernhard et al., 2022b). Based on the assumptions that (1) ground ice starts at the depth where the active layer ends and (2) SOC is not present in any form in the massive ground ice, we sampled the negative elevation change within each RTS label to calculate the SOC mobilisation per RTS

$$\begin{array}{}\text{(3)}& \begin{array}{rl}{\mathrm{SOC}}_{\mathrm{RTS}}& =\sum ^{n\left(\mathit{\delta }A\right)}\sum _{d=\mathrm{1}}^{n=\mathrm{10}}\left[\mathrm{SOC}\right(d_{<\mathrm{ALT}})\\ & +\mathrm{SOC}\left(d_{>\mathrm{ALT}}\right)(\mathrm{1}-\text{GI}/\mathrm{100})]\mathrm{100}{\mathrm{m}}^{\mathrm{2}}\end{array}\end{array}$$

with the number of pixels n, the RTS ablation area δA [m], the depth of the active layer ALT (>0 m), ground ice content GI [%] and soil depth d ($\mathrm{0}<d<\mathrm{10}$ m). If the depth d is greater than the ALT, then only the part of the eroded material that is not massive ground ice is added to the total SOC mobilisation. Similarly, we estimate the volume of RTS-induced ground ice loss across the depth layers by scaling the eroded material by $\mathrm{1}-\text{GI}/\mathrm{100}$. To illustrate, if we sample at an exemplary RTS location an ALT = 1.8 m and GI${}_{\mathrm{2}-\mathrm{3}\mathrm{m}}=\mathrm{31}$ % we round the ALT to 2 m. For 0–2 m depth, we compute the SOC mobilisation without scaling for ground ice since we assume no presence of massive ice in the active layer. For the depth layer 2–3 m we reduce the SOC mobilisation by 31 %.

We report the total estimates of the volume of eroded material, the mobilisation of SOC, and the loss of ground ice as the sum over all RTS in the study region and throughout the study period. By dividing SOC mobilisation by the number of years between T1 and T2, we estimate a yearly SOC mobilisation rate. However, the values might be partially biased since we cannot distinguish between RTS that have been active for the entire study period and those that may have only been active for a shorter period of time (Bernhard et al., 2022b). We normalise our results by the size of the study area and the number of RTS to ensure comparability to the results of previous studies.

2.2.3 RTS material erosion error assessment and validation

Since we aim to integrate two datasets with different spatial resolutions (DEM raster: 10 m, optical RTS inventory vectors: 3 m), we have to account for the ambiguity produced. Negative elevation change pixels might only partially intersect the RTS delineation instead of being entirely contained. We chose to rasterise the optical RTS delineations into two parts: an upper-bound RTS area δA⁺ that includes all intersected pixels, and a lower-bound area δA⁻ with only fully contained pixels. Together with the estimated elevation error σ_h of the DEM we compute the upper and lower volume change bounds

$$\begin{array}{}\text{(4)}& \mathit{\delta }{V}^{+}=\sum _{n\left(\mathit{\delta }{A}^{+}\right)}(\mathit{\delta }{h}_{<\mathrm{0}\mathrm{m}}-{\mathit{\sigma }}_{\mathrm{h}})\cdot \mathrm{100}{\mathrm{m}}^{\mathrm{2}}\end{array}$$

and

$$\begin{array}{}\text{(5)}& \mathit{\delta }{V}^{-}=\sum _{n\left(\mathit{\delta }{A}^{-}\right)}(\mathit{\delta }{h}_{<\mathrm{0}\mathrm{m}}+{\mathit{\sigma }}_{\mathrm{h}})\cdot \mathrm{100}{\mathrm{m}}^{\mathrm{2}}\end{array}$$

that indicate the uncertainty induced by differences in spatial resolution, boundary mismatch, and vertical DEM error. Together with the uncertainty estimates of the SOC stocks, we propagate the error bounds to all reported quantities. We assume the RTS area A_Xia is error-free, despite biases and subjective influences in both automated and manual RTS segmentation, which are difficult to measure (Nitze et al., 2024; Maier et al., 2025).

In addition to reporting data uncertainties, we performed a two-fold validation at five validation sites (Fig. 1a): (1) We manually delineated the RTS ablation area in the elevation change maps, where a distinct pattern of negative elevation change is visible to the human eye. At all validation sites, the elevation change maps consisted of TanDEM-X DEMs from 2011 and 2019. We statistically compared the delineations of the RTS inventory (Xia et al., 2024) from the same year (2019), the year before (2018) and the year after (2020) to the manually delineated ablation zones to investigate the agreement between the datasets in terms of RTS quantities, (ablation) area, and material erosion volume (Fig. 2d). (2) At the Beiluhe River Basin site, we compared the TanDEM-X-derived elevation change with very high-resolution (VHR) photogrammetric DEMs from an in-situ drone campaign covering in total six RTSs in August 2020. A DJI P4 Multispectral was used to obtain the multispectral drone images. The resulting DEMs have a spatial resolution of <1 m and a georeferencing accuracy of 0.2 m RMSE. Since no VHR DEM was available for T1, we could not perform DEM differencing and directly validate our volume change estimates. Based on the hillshade VHR DEMs, we manually delineated the approximate location of the headwall with the help of transect profiles (Fig. 2d). We defined small buffer zones (∼5 m) and randomly distributed points (n=100 per RTS) on both sides of the headwalls that represent the elevation of stable ground h_stable and the RTS slump floor h_RTS, respectively. We computed the average headwall height h_VHR per RTS as the median difference between h_stable and h_RTS. The monitored RTS were relatively small (<10⁴ m²) and shallow (h_VHR<4 m). Defining a headwall position and applying the same methodology with TanDEM-X DEMs is not feasible due to the coarser resolution. Therefore, we chose to compare the estimated headwall heights based on the VHR DEMs with the maximum negative elevation change δh_max that we estimated based on TanDEM-X elevation change maps, assuming that the largest height loss aligns with the largest material ablation and is located close to the headwall.

3 Results

3.1 RTS mass wasting on the QTP between 2011 and 2020

Out of 3613 RTS delineations of Xia et al. (2024) that we matched with the generated elevation change maps, we excluded 2 % (71 RTS) from further analysis due to low SAR coherence or SAR layover/shadow regions. For these discarded RTS, reliable erosion volume estimates were not possible to achieve due to the low confidence in DEM quality.

Figure 3 Total RTS mass wasting on the QTP between 2011–2020. All data are aggregated in 50×50 km-tiles for visualisation purposes. For all mass wasting quantities including (a) material erosion volume, (b) total GI loss, and (c) total SOC mobilisation, the central QTP shows the highest values.

We estimated a total volume change of eroded material δV of 5.02${}_{\mathrm{0.75}}^{\mathrm{25.35}}\times {\mathrm{10}}^{\mathrm{7}}$ m³ induced by RTS activity on the QTP between 2011 and 2020 (Fig. 3a). Approximately half of the volume change originates from 0–1 m, 28 % from 1–2 m, and 13 % from 2–3 m depth. On average, 65 % of the entire delineated areas A_Xia of the RTS inventory were actively eroding between 2011 and 2020. The median ablation area δA per RTS was 5200 m² compared to the full delineation A_Xia (including the accumulation and inactive parts of an RTS) of 8000 m². We estimated a median volume loss of 6534±2284 m³. When fitting a linear model to the log-transformed area δA and material erosion volume δV based on (1) OLS and (2) ODR, we found power-law relationships for the area-volume scaling of RTS on the QTP (Fig. S1a) of

$$\begin{array}{}\text{(6)}& \begin{array}{rl}& \mathit{\delta }{V}_{\mathrm{ODR}}=(\mathrm{0.09}\pm \mathrm{0.01})\cdot \mathit{\delta }{A}^{\mathrm{1.30}\pm \mathrm{0.01}}\\ & \text{with }{R}^{\mathrm{2}}=\mathrm{0.87}(p<\mathrm{0.001})\end{array}\end{array}$$

and

$$\begin{array}{}\text{(7)}& \begin{array}{rl}& \mathit{\delta }{V}_{\mathrm{OLS}}=(\mathrm{0.22}\pm \mathrm{0.01})\cdot \mathit{\delta }{A}^{\mathrm{1.20}\pm \mathrm{0.01}}\\ & \text{with }{R}^{\mathrm{2}}=\mathrm{0.87}(p<\mathrm{0.001}).\end{array}\end{array}$$

For ODR, we obtain the same scaling coefficient (α_ODR=1.30, c=0.05, Fig. S3a) when we use the entire area of the RTS delineations A_Xia instead of solely the ablation area δA yet a lower scaling coefficient for computations based on OLS (α_OLS=1.11, c=0.29, Fig. S3a). However, the fit is slightly noisier (R²=0.75 for ODR and R²=0.77 for OLS). An α value between 1.11 and 1.30 indicates that RTS on the QTP followed a relationship between a growing scar zone with constant depth (α=1.0) and growth with a constant width-depth ratio (α=1.5) during the last decade and fall in the range of soil landslides (1.1–1.4) based on the investigated scaling relations of landslides in temperate climate regions (Jaboyedoff et al., 2020; van der Sluijs et al., 2023).

Figure 4 RTS mass wasting quantities in the QTP subregions between 2011 and 2020: The first row displays box plots for values associated with individual RTS, the second row bar charts for the total quantities across the subregions, and the third row additional data that vary between the columns. For all plots, the computed uncertainty is reported with error bars. (a) RTS material erosion and area-volume scaling coefficients for the subregions with R² values between 0.77/0.79 (p<0.001, West) and 0.89/0.90 (p<0.001, Central) computed based on ODR (black dots)/OLS (purple dots). (b) RTS GI loss and the distribution of GI content (Zou et al., 2024). (c) RTS-induced SOC mobilisation and the distribution of SOC stocks 0–3 m (Wang et al., 2021). Number of RTS per subregion: West = 170, West-Central = 523, Central = 2688, East = 76, Northeast = 140.

We estimated that ${\mathrm{3.58}}_{\mathrm{0.28}}^{\mathrm{28.20}}\times {\mathrm{10}}^{\mathrm{6}}$ m³ of massive ground ice has been lost on the entire QTP during the last decade (Fig. 3b). Approximately 64 % of the thawed ground ice was located in the first metre under the active layer (2–3 m), 32 % between 3 and 5 m and the remaining 4 % below 5 m depth. Based on the SOC stocks dataset (Wang et al., 2021), we calculated a total SOC mobilisation of ${\mathrm{2.78}}_{\mathrm{0.11}}^{\mathrm{7.98}}\times {\mathrm{10}}^{\mathrm{8}}$ kg C (Fig. 3c) and annual SOC mobilisation rates of 0.35×10⁸ kg C a⁻¹. The first metre of soil contributed approximately 76 %, the second 14 %, and the third 8 % to the total SOC mobilisation.

Xia et al. (2024) found 75 % of RTS activity in the central QTP, including the highest area expansion rates. Figure 4 shows a similar pattern of 78 % of material erosion volume, 89 % of ground ice loss and 81 % of SOC mobilisation attributed to the central QTP. We found that even though the central Plateau exhibited the highest absolute amounts due to the highest density of RTS, individual RTS in the northeastern QTP eroded on average more material (Fig. 4a). The scaling coefficients range from α_ODR=1.27–1.34 (R²=0.77–0.89, p<0.001) and α_OLS=1.11–1.23 (R²=0.79–0.90, p<0.001) in the West to East subregions to ${\mathit{\alpha }}_{\mathrm{ODR}}=\mathrm{1.47}\pm \mathrm{0.05}$ (R²=0.87, p<0.001) and ${\mathit{\alpha }}_{\mathrm{OLS}}=\mathrm{1.34}\pm \mathrm{0.04}$ (R²=0.87, p<0.001) in the northeast QTP. Based on Zou et al. (2024), the GI content is highest in the central subregions (median of 32.8 % in central and 32.6 % in the west-central QTP), where we also found the highest total and average ground ice loss (Fig. 4b). The lowest GI content was present in the northeast (median of 0 %, mean of 10.7 %). The total amount of SOC mobilised from RTS activity for the QTP subregions showed a pattern similar to material erosion and ground ice loss: the central QTP dominated all RTS mass wasting quantities through the largest number of RTS present. The average SOC mobilisation per RTS increased from the west to the northeast QTP (Fig. 4c). Details can be found in the Supplement (Table S1a, b and Fig. S2).

3.2 Assessment of InSAR DEM- and optical-based monitoring of RTS mass wasting on the QTP

For the validation sites, we compared the delineations (2018–2020) of the RTS inventory (Xia et al., 2024) to the ablation zones we manually identified on the elevation change maps based on TanDEM-X DEMs (2011–2019). Out of a total of 445 RTS in the RTS inventory, we identified 290 in the lower-resolution elevation change maps, which accounts for an F1 score of 0.63. 17 RTS were missed in the RTS inventory, while 155 RTS present in the RTS inventory were not distinguishable from background noise in the elevation change maps. Most of the validation sites showed a good agreement between the two datasets. However, in Western Kunlun (A in Fig. 1a) only 55 % of the RTS in the RTS inventory were detectable in the elevation change maps. More details on individual results at the test sites can be found in Table S2a and b.

Figure 5 Compatibility and accuracy assessment aggregated for all validation sites: (a) RTS at the Beiluhe River Basin site in an high-resolution optical image from 2018 (ESRI Satellite Basemap) and on a TanDEM-X elevation change map (2011–2019) with delineations (2018–2020) from the RTS inventory (Xia et al., 2024) (solid lines). The manually delineated ablation area (= negative elevation change) is visualised by a dashed line. The RTS has grown over the course of the three years and its headwall extended upslope. The delineations from the RTS inventory based on optical images and disturbances of the vegetation cover include not only ablation zones but also material accumulation further downslope (= positive elevation change). (b) Sum of ablation and accumulation area based on the elevation loss/gain pixels for the delineations of the RTS inventory. Delineations of the RTS inventory tend to cover a larger area than actual ablation area distinguishable on the DEM. (c) Distribution of the RTS ablation area. Due to the higher resolution of the optical images compared to the TanDEM-X DEMs, smaller RTS can be distinguished from the image background. Only small differences can be observed between the years. (d) Sum of material erosion volume based on the delineations of the RTS inventory and the negative elevation change. The volume computed from the 2019 delineation was closest to the actual erosion volume while 2018 under- and 2020 overestimates the actual erosion volume. (e) For six RTS at the Beiluhe River Basin site (D in Fig. 1a), we compared the maximum elevation loss δh_max of the TanDEM-X elevation change within the 2020 delineation of the RTS inventory to the average headwall height derived from drone-based single-time-step VHR DEMs (summer 2020). We assume that the maximum elevation loss can be used as an approximation of the headwall height of an RTS.

Visually, the delineations from the RTS inventory fit well with the spatial patterns of the TanDEM-X elevation change map, with most of the ablation zone being covered. Due to the lower spatial resolution of the DEM, the differences between the delineations from the RTS inventory of 2018, 2019, and 2020 are rather small in this example (Fig. 5b). The total RTS area A_Xia consists of the material erosion zone (= negative elevation change) which was growing between 2018 and 2020, and the material accumulation zone (= positive elevation change) which stayed relatively constant across time. For all analysed years, the total RTS area was distinctly larger than the sum of the manually delineated erosion areas. When only considering the ablation zone within the RTS delineation, the difference was smaller (Fig. 5b). The actual average RTS ablation area based on the manual delineation was 1.29×10⁴ m², which is relatively close to the RTS inventory: 0.98×10⁴ m² in 2018, 1.13×10⁴ m² in 2019, and 1.35×10⁴ m² in 2020, accounting for 62 % to 64 % of the total RTS area A_Xia. However, in the high-resolution PlanetScope images used as the basis for the RTS inventory, considerably smaller RTS could be identified compared to the TanDEM-X DEMs. We found that RTS monitoring using TanDEM-X can only identify RTS with sizes larger than 10³ m² (Fig. 5d). However, in our analysis of the entire QTP, only 6 % of 3613 RTS had ablation areas below this value. The total volume of material erosion computed with the RTS inventory differed only minimally from the volume calculated based on the manual delineations (Fig. 5e). We estimated a total volume δV of 8.47×10⁶ m³, 9.98×10⁶ m³, and 10.97×10⁶ m³ for the RTS inventory delineations of 2018, 2019, and 2020, respectively, compared to an erosion volume of 9.94×10⁶ m³ defined by the manual delineations. The uncertainty of the estimated material erosion volume is larger for the results based on the delineations of the RTS inventory compared to the results based on manually delineated ablation areas (Table S2b).

We used VHR DEMs to further validate our results. Figure 5f shows an adequate fit between the maximum elevation changes δh_max computed with TanDEM-X DEMs and the delineations of the RTS inventory and the average headwall height calculated from the VHR DEMs at six RTS locations. Details can be found in Table S3. However, the small sample size does not allow meaningful statistical analysis and, therefore, only allows for a qualitative comparison.

4 Discussion

4.1 RTS activity, material erosion, and area-volume scaling across permafrost regions

Previous research has investigated material erosion volumes (Lantuit and Pollard, 2005; Kokelj et al., 2015; Günther et al., 2015) and allometric scaling relationships for thaw-driven mass wasting primarily in regions <10⁴ km² in the Arctic (Kokelj et al., 2021; van der Sluijs et al., 2023) and recently on the QTP (Ma et al., 2025). At a larger spatial scale, comparable to our study, Bernhard et al. (2022a) estimated annual material erosion rates for ten sites across Canada, Alaska, and Siberia, covering a total area of 220 000 km², using DEMs derived from TanDEM-X observations. Similarly, Dai et al. (2025) used ArcticDEM time series to assess volumetrics and area-volume scaling across other Arctic sites. On the QTP, we found RTS material erosion rates of 6.36${}_{\mathrm{0.95}}^{\mathrm{32.1}}\times {\mathrm{10}}^{\mathrm{6}}$ m³ a⁻¹, whereas Bernhard et al. (2022a) estimated rates ranging from $<\mathrm{0.5}\times {\mathrm{10}}^{\mathrm{6}}$ m³ a⁻¹ in Alaska (Fig. 6, location 5) and the Siberian Arctic (6–8) to 7.16×10⁶ m³ a⁻¹ in the Canadian Arctic (1–4) between 2010 and 2017. With 5.2 RTS per 100 km² between 2010 and 2020, the QTP has an RTS density approximately half of that observed between 2010 and 2017 on Banks Island (2) (Bernhard et al., 2022a), yet more than double that of all other Arctic sites investigated. Over the past decade, RTS on the QTP displaced 101.8 ${\mathrm{m}}^{\mathrm{3}}{\mathrm{a}}^{-\mathrm{1}}{\mathrm{km}}^{-\mathrm{2}}$ material, which is roughly 4, 6, and 10 times less than the Canadian hillslope RTS hotspots, the Peel Plateau (4), Ellesmere Island (1), and Banks Island (2), respectively (Lewkowicz and Way, 2019). The remaining sites investigated by Bernhard et al. (2022a) in the Canadian (3) and Siberian (6–8) Arctic are characterised by smaller and shallower RTS located on lake shores in relatively flat terrain and exhibited less than half the volume change rates observed on the QTP (Fig 6b).

Figure 6 Comparison between RTS on the QTP and in Arctic permafrost regions with respect to material erosion and area-volume scaling. (a) Map of the Arctic sites investigated by previous studies: Kokelj et al. (2021) and van der Sluijs et al. (2023) estimated erosion volume and area-volume scaling (OLS) based on high-resolution airborne DEMs and pre-disturbance reconstruction in northwestern Canada. Bernhard et al. (2022a) studied RTS area-volume relations based on TanDEM-X DEMs (ODR) at eight sites in North America and Siberia between 2010 and 2016. (b) RTS density and material erosion volume per unit area for all sites reported in Bernhard et al. (2022a). RTS density and erosion volume are consistently higher on the QTP compared to Siberia and in a similar magnitude as the North American sites. (c) The area-volume scaling coefficients (ODR) reported by Bernhard et al. (2022a) are on a similar magnitude as α_ODR of the QTP. The Siberian sites have generally lower α-values while Banks Island (2) and Peel Plateau (4) are closest to the QTP's α_ODR. Both Canadian sites (2, 4) are dominated by hillslope RTS compared to a prevalence of lakeshore RTS in Siberia (6–8). The coefficients based on OLS reported by Kokelj et al. (2021) and van der Sluijs et al. (2023) are distinctly higher than α_OLS on the QTP.

Area-volume or allometric scaling relationships describe volumetric enlargement characteristics of RTS and their potential drivers such as that higher α coefficients indicate larger headwalls and concavity depth per unit area growth (van der Sluijs et al., 2023). Several studies investigated the power-law relationship between RTS area and volume using different methodologies and datasets, complicating direct transferability between results. Similar to our approach, Bernhard et al. (2022a) and Dai et al. (2025) calculated elevation change over a time period T1–T2, whereas Kokelj et al. (2021) and van der Sluijs et al. (2023) derived erosion volumes by subtracting a simulated pre-disturbance DEM (T1) from the high-resolution DEM of the disturbed state (T2) using masking and re-interpolation techniques. Ma et al. (2025) followed a similar approach yet used 30 m-resolution mosaic DEMs as the pre-disturbance state (T1). The studies also used different model fitting approaches: OLS (Kokelj et al., 2021; van der Sluijs et al., 2023; Dai et al., 2025) and ODR (Bernhard et al., 2022a). Ma et al. (2025) did not explicitly report the model choice. We therefore assume that the authors used the more common OLS. Bernhard et al. (2022a) estimated an overall Arctic ${\mathit{\alpha }}_{\mathrm{ODR}}=\mathrm{1.15}\pm \mathrm{0.01}$ (n_RTS=1854) while Dai et al. (2025) found with ${\mathit{\alpha }}_{\mathrm{OLS}}=\mathrm{1.30}\pm \mathrm{0.01}$ (n_RTS=2747) a noticeably higher coefficient. Kokelj et al. (2021) and van der Sluijs et al. (2023) estimated α_OLS-values of 1.41 (n_RTS=71) and 1.36±0.01 (n_RTS=1522) in the low Canadian Arctic, respectively (Fig. 6c, 3 + 4). Bernhard et al. (2022a) calculated the power-law scaling separately and found ${\mathit{\alpha }}_{\mathrm{ODR}}=\mathrm{1.26}\pm \mathrm{0.02}$ (n_RTS=438) for the Peel Plateau (3) and ${\mathit{\alpha }}_{\mathrm{ODR}}=\mathrm{1.16}\pm \mathrm{0.03}$ (n_RTS=212) for the Tuktoyaktuk Coastlands (4). This range of scaling coefficients for similar regions as well as our results for QTP (Fig. 4a and Fig. S1) highlights the challenge in comparing scaling studies based on different methodologies and datasets. Moreover, minor differences in the scaling coefficient have strong impacts on the scaling: A difference of, for example, 0.1 in α leads to a doubling of the volume estimation (van der Sluijs et al., 2023).

On the QTP, more than half of the RTS in the inventory of Xia et al. (2024) initiated in 2016 (Fig. S2a) while Luo et al. (2022) found that >80 % of the RTS in the central QTP formed during extremely warm summers in 2010 and 2016. This indicates that our scaling results are based on planimetric areas and volumes representing rather the entire RTS landform and its lifecycle, similar to the methodology used by Kokelj et al. (2021), van der Sluijs et al. (2023) and Ma et al. (2025). Compared to Arctic permafrost regions, our results (${\mathit{\alpha }}_{\mathrm{ODR}}=\mathrm{1.30}\pm \mathrm{0.01}$, ${\mathit{\alpha }}_{\mathrm{OLS}}=\mathrm{1.20}\pm \mathrm{0.01}$, n_RTS=3043) are higher than the scaling coefficients reported by Bernhard et al. (2022a) yet lower than the estimates of Dai et al. (2025), van der Sluijs et al. (2023), and Kokelj et al. (2021). When comparing material erosion and area-volume scaling in the QTP subregions, we see almost as pronounced differences as between the QTP and the Arctic (Fig. 4a). The central QTP dominates the number of RTS and the associated total material erosion (Xia et al., 2024). However, the northeast QTP has a substantially higher scaling coefficient, indicating larger headwalls and hence, potentially more efficient material mobilisation at an individual RTS level. The mountainous northeastern plateau potentially favours the development of relatively larger and deeper RTS (see van der Sluijs et al., 2023, for coefficient interpretation, Figs. S1f and S3f). However, the driving factors underlying the observed differences across permafrost (sub-)regions are likely multifaceted, including variations in RTS longevity and lifecycle stage, GI content, vegetation and soil properties, and proximity to water bodies. Detailed investigations into these factors remain the subject of prospective research. RTS on the QTP are reported to be generally smaller and more shallow than at Arctic hotspots (Liu et al., 2024). It is therefore possible that the coarse resolution of the TanDEM-X DEM might not correctly capture the area and volume change for these small areas, skewing the scaling models. However, using high-resolution stereo-optical DEMs, Ma et al. (2025) found a scaling coefficient similar to our α_OLS of 1.20±0.01 (n_RTS=1429). Area-volume scaling also depends on how the actively eroding RTS scar zone is delineated. While there are commonly accepted definitions of RTS geomorphology (CPA et al., 2024; Harris, 1988), different data sources show different aspects of RTS activity, and even experts working on similar data can differ strongly in their delineations based on their own ontological understandings of what constitutes RTS (Nesterova et al., 2024; Nitze et al., 2024; Maier et al., 2025). More research is needed in this regard, and our novel dataset offers a critical resource for future investigations into the mechanisms driving RTS material erosion on the QTP. Our results are consistent with other studies showing a notable increase in RTS (mass wasting) activity across the QTP, particularly over the past decade (Ma et al., 2025; Xia et al., 2022; Luo et al., 2022) with RTS sizes being generally smaller and headwall retreat rates lower than those found in other Arctic regions (Yi et al., 2025; Luo et al., 2022; Lewkowicz and Way, 2019; Runge et al., 2022; Nesterova et al., 2024; Huang et al., 2021). The magnitude of newly formed RTS potentially offsets the relatively low concavity depths, so that the total mass-wasting activity and material erosion volume during the last decade show a magnitude comparable to the thermokarst landscapes in the high Arctic.

4.2 Magnitude of RTS SOC mobilisation and ground ice loss

Jiao et al. (2022) investigated one large RTS in the Beiluhe River Basin between 2021 and 2022 and found 1.9 m a⁻¹ vertical deformation at the headwall and a total volume change of 1.41×10³ m³ a⁻¹. The active layer at the RTS location was 1.95 m with a ground ice layer between 2.2 and 3.5 m depth and an ice content of 68 % to 88 % at depths of 2.2 to 4 m obtained at a borehole near the RTS. Compared to our results in the central QTP (Fig. 4, Table S1a), the investigated RTS has a typical headwall height and material erosion volume ($\stackrel{\mathrm{‾}}{\mathit{\delta }h}=\mathrm{1.45}$ m, δh_max=2.35 m, $\stackrel{\mathrm{‾}}{\mathit{\delta }V}=\mathrm{1.97}\times {\mathrm{10}}^{\mathrm{3}}$ m³ a⁻¹). The measured GI content near the RTS is substantially higher than the average GI content we used in our estimations ($\stackrel{\mathrm{‾}}{\mathrm{GI}}=\mathrm{32.4}$ %). The GI content is strongly related to the landscape's geomorphology. The majority of GI on the QTP is found on gentle shaded slopes at elevations between 4400 and 5100 m (Fan et al., 2023) typical for central subregions and most RTS locations (Table S1, Fig. S2b–d). However, some studies also found GI content higher than 80 % in the northeast QTP (Wang et al., 2018; Fan et al., 2023) and a mean GI content of ∼ 16 % in the Beiluhe River Basin (Lin et al., 2020). RTS only form in locations where massive ground ice is present in depths that can be exposed by, for example, active layer detachments (Nesterova et al., 2024). We might underestimate local ground ice conditions since large-scale datasets that estimate soil properties in permafrost regions are based on limited observations and coarse spatial resolution that typically do not well represent fine-scale soil conditions (Hugelius et al., 2014; Mishra et al., 2021; Wang et al., 2021; Zou et al., 2024). With increased soil depth, estimates of soil properties become even more scarce (Chen et al., 2024b; Ding et al., 2019). Since we observed negative elevation changes of more than 3 m in RTS on the QTP, we used a simple exponential model to extrapolate SOC stocks to deeper soil layers. Probably, this model is too simple and does not capture the spatial variability of soil conditions in the complex permafrost landscape (Bernhard et al., 2022b). However, only 2 % of the total estimated mobilised SOC came from soils below 3 m depth.

The processes and fate of the mobilised SOC remain highly uncertain as multiple complex ecosystem interactions and hydrothermal processes are involved (Yi et al., 2025). Parts of mobilised SOC remain potentially on the slump floor and are available for microbial decomposition and release as greenhouse gases (Wang et al., 2024). Other parts, together with the thawed material, are deposited and stabilised at the slump floor (Thomas et al., 2023; Liu et al., 2021, 2018; Mu et al., 2017) or laterally transported downslope into adjacent river and lake systems undergoing complex water chemistry processes such as dissolution or sedimentation (Lewkowicz and Way, 2019). We only quantified the magnitude of SOC mobilised by RTS activity and its spatial pattern which closely follows the spatial trend of existing SOC stocks on the QTP (Wang et al., 2021). Exploring the complex pathways of mobilised SOC is beyond the scope of this study. However, based on the results of our study, we recommend further research into the fate of the mobilised SOC. Figure 5b clearly shows areas of material ablation, but also regions of material deposition on the slump floor are visible. In future investigations, insights could potentially be generated into how much relocated material and SOC are deposited within an RTS and how much is transported laterally into hydrological networks. In addition, the balance between areas of negative and positive elevation change might be another avenue to improve spatially explicit GI information for the QTP.

To our knowledge, no study has quantified the SOC mobilisation based on empirical RTS erosion material volumetrics for the entire QTP. Ma et al. (2025) modelled a total annual SOC loss of 4.12×10⁷ kg C a⁻¹ (95 % CI: 3.06×10⁷ kg C a⁻¹–5.12×10⁷ kg C a⁻¹) from RTS mass wasting between 1989 and 2022 based on optical RTS inventories and allometric scaling relations, which is of a similar magnitude as our results (${\mathrm{3.53}}_{\mathrm{0.15}}^{\mathrm{10.13}}\times {\mathrm{10}}^{\mathrm{7}}$ kg C a⁻¹). Bernhard et al. (2022b), similar to this study, processed multi-temporal stacks of TanDEM-X-derived DEMs for the time periods 2010–2016 (n_RTS=76) and 2017–2020 (n_RTS=1404) to investigate the impact of a summer heat wave on RTS mass wasting and SOC mobilisation on the Taymyr Peninsula in Siberia. Individual RTS on the Taymyr Peninsula mobilised on average a similar amount of material, yet 14 times more SOC compared to the RTS on the QTP prior to the heatwave. The heatwave triggered a surge in RTS activity, including a 17-fold increase in new initiations, a 2.3-fold increase in the average volume change rate, and a 28-fold increase in SOC mobilisation, resulting in quantities far exceeding those observed on the QTP.

Due to the increase in greening and wetting of the QTP, recent studies estimate that the QTP is a large carbon sink of 344×10⁸ kg C a⁻¹ (or 34.4 Tg C a⁻¹) (Chen et al., 2024a; Wang et al., 2023). However, several studies that conducted soil sampling within disturbed permafrost areas for several consecutive years found that up to one-third of the surface SOC content (< 40 cm) has the potential to be lost due to rapid permafrost thaw (Mu et al., 2017; Liu et al., 2018; Wu et al., 2018). Furthermore, vegetation restoration appears to be slow on the QTP, and exposed areas disturbed by thermokarst can remain bare for decades (Liu et al., 2018; Mu et al., 2017; Li et al., 2025). Even though our results indicate that RTS-induced SOC mobilisation only accounted for ∼0.1 % of the total QTP carbon budget during the last decade and is proportionately insignificant to many Arctic regions, the sharply rising number of RTS on the QTP and the majority of SOC mobilisation in the uppermost, carbon-rich layers of soil could become more relevant for regional carbon budgets in a warming future climate.

4.3 Limitations and Potential of a Multimodal Data Approach for RTS Mass Wasting Monitoring

Compared to the Arctic where the ArcticDEM strip data (Porter et al., 2018) offers an open source and high-resolution multitemporal DEM source based on stereo-optical satellite images (∼ 3 m) that is suitable for RTS monitoring (Dai et al., 2024, 2025; Nitze et al., 2021; Yang et al., 2023), on the QTP no similar high-resolution temporally resolved DEM exists - except for the data from the TanDEM-X mission. Bistatic TanDEM-X observations uniquely enable DEM generation on a global scale. Even in regions with adverse geographic or climatic conditions for satellite remote sensing, for example, a high percentage of cloud cover or long periods of snow cover, bistatic radar observations can be used to produce high-quality DEMs with acceptable spatial resolution for RTS monitoring (Krieger et al., 2007; Bojarski et al., 2021; Bernhard et al., 2022a; Maier et al., 2025). However, the temporal resolution of the TanDEM-X observations is not equal for all regions of the world. We can assume that RTS activity only occurs in the warm summer months and a phase of stability occurs between October and April due to low temperatures (Chen et al., 2015; Che et al., 2008; Ma et al., 2023). Due to limited data availability on the QTP, we had to aggregate observations from several years and seasons to ensure sufficient coverage. For T1, we used observations between mid-winter 2011 and end of winter 2012 accounting for exactly one summer of RTS activity, while for T2 we had to accept a larger time span (end of summer 2017 to mid-winter 2020) accumulating the RTS activity of two summer seasons, which potentially increased imprecision in the mass wasting estimates. Xia et al. (2024) found that the highest RTS activity occurred before 2020. Between 2020 and 2022, only 59 new RTS were detected in the PlanetScope images. This could indicate that we captured most of the RTS activity that has occurred on the QTP during the last decade. For volume estimation based on delineations of the RTS inventory, stable areas containing no elevation change but background noise are likely included, since most delineations based on optical images are broader than the actual active erosion area. Although this minimally affects the total volume change due to the low magnitude in negative elevation change, it adds additional random errors, thus contributing to the overall uncertainty budget.

To monitor the dynamic lifecycles of these complex thermokarst features and understand their drivers and future development, yearly records of RTS-induced material erosion volumes are highly desirable (Nesterova et al., 2024; Kokelj et al., 2021). The general scarcity of suitable DEMs, as well as the temporal limitations, makes monitoring of RTS material erosion and impact on carbon cycles solely based on DEM data challenging. Allometric or area-volume scaling relations, as presented and discussed in this study, are typically used to investigate landscape evolution or RTS activity change over time, but can also be used to enable the investigation of not only planimetric area expansion from optical RTS inventories but also transform the area change into volume change of eroded material and expand the analysis to yearly or even seasonal temporal scales. We derived a statistically robust area-volume scaling relation for the entire QTP region, similar to regions elsewhere (Bernhard et al., 2022a; van der Sluijs et al., 2023; Dai et al., 2025). However, a small difference in the scaling coefficient α has a large impact on the resulting volume change δV. We see substantial differences in the scaling relations between the different QTP subregions (Figs. S1 and S3). When comparing our scaling results with the small number of existing studies on allometric scaling of RTS in the Arctic, we saw distinct differences between our results and studies investigating similar regions in the Arctic (Bernhard et al., 2022a; van der Sluijs et al., 2023; Kokelj et al., 2021) yet clear similarities to the only existing study on the QTP (Ma et al., 2025). Applying area-volume scaling especially to multimodal RTS datasets should be done carefully to obtain regional estimates on material erosion volume and mass wasting derivatives such as ground ice loss or SOC mobilisation. For subregional scale or even feature level, this approach has its clear limitations due to the heterogeneity of RTS across at both local and broad scales. However, large-scale and even pan-Arctic RTS inventories based on optical satellite images become more available and may allow for a similar approach of finding scaling relations for Arctic thermokarst regions if delineations coincide well ontologically, spatially and temporally with measured elevation change (Yang et al., 2025; Nitze et al., 2025). Here, we present a scaling relationship for the QTP based on commonly applied OLS regression

$$\begin{array}{}\text{(8)}& \begin{array}{rl}& \mathit{\delta }{V}_{\mathrm{OLS}}=(\mathrm{0.29}\pm \mathrm{0.01})\cdot {A_{\mathrm{Xia}}}^{\mathrm{1.11}\pm \mathrm{0.01}}\\ & \text{with }{R}^{\mathrm{2}}=\mathrm{0.77}(p<\mathrm{0.001})\end{array}\end{array}$$

that may readily transform planimetric RTS area derived from vegetation disturbance on optical remote sensing images into volume estimates at scale. We showed that RTS monitoring on elevation change maps based on DEMs with 10 m resolution omits approximately 35 % of the present features compared to monitoring RTS with 3 m multispectral images. However, the difference in estimated material erosion volume from the two datasets is <1 % since most of the missed RTS are small and shallow (Fig. 5c). Open-source and multi-temporal images, such as from ESA's Sentinel-2 satellites, with a similar spatial resolution as TanDEM-X DEMs could have great potential to continuously monitor RTS activity to answer questions about volumetrics, permafrost thaw impacts on hydrological systems, and carbon cycles based on reliable area-volume scaling laws.

5 Conclusions

RTS landforms are typically complex and highly dynamic, often remaining active for several years before stabilising and, in some cases, reinitiating. To adequately capture their temporal evolution and interactions with changing climate, high-temporal-resolution remote sensing data are essential. By combining modelled soil property datasets with multimodal remote sensing data, we estimated that RTS activity relocated 5.02${}_{\mathrm{0.75}}^{\mathrm{25.35}}\times {\mathrm{10}}^{\mathrm{7}}$ m³ formerly frozen material, contributed to a loss of 3.58${}_{\mathrm{0.28}}^{\mathrm{28.20}}\times {\mathrm{10}}^{\mathrm{6}}$ m³ ground ice and mobilised 2.78${}_{\mathrm{0.11}}^{\mathrm{7.98}}\times {\mathrm{10}}^{\mathrm{8}}$ kg C SOC on the QTP between 2011 and 2020. Interregional comparisons of RTS dynamics are challenging due to varying spatial, temporal, and methodological factors. However, RTS on the QTP exhibit mass-wasting dynamics comparable to Arctic RTS hotspots. Despite their comparably recent initiation and smaller size, erosion and SOC mobilisation on the QTP in the past decade surpassed some regions in the Siberian Arctic but remained up to 10 times lower than well-known thermokarst regions in the high Canadian Arctic. Although RTS-induced carbon mobilisation only accounts for approximately 0.1 % of the QTP's carbon budget, the acceleration of RTS activity since the beginning of the last decade can contribute to the anticipated region's shift from a carbon sink to a source. By integrating remote sensing data with varying spatial and temporal resolutions and different information layers, we demonstrated that erosion volumes can be accurately estimated, even when the delineations of the RTS erosion-affected area vary. We found a reliable power-law scaling based on commonly applied OLS regression between the computed material erosion volume change and the (ablation) area in the optical RTS inventory ($\mathit{\delta }V=\mathrm{0.05}\cdot \mathit{\delta }{A}^{\mathrm{1.20}\pm \mathrm{0.01}}$, Fig. S1a, and $\mathit{\delta }V=\mathrm{0.29}\cdot {A_{\mathrm{Xia}}}^{\mathrm{1.11}\pm \mathrm{0.01}}$, Fig. S3a) that enables future research to transform the planimetric area of RTS delineations into regional estimates of erosion volume and constrain RTS-induced SOC mobilisation on the QTP. Improved estimates and allometric relationships will help close the knowledge gap in understanding the impact of permafrost thaw on the permafrost carbon cycle for the QTP and globally.