A Comparative Study of Active Rock Glaciers Mapped from Geomorphic- and Kinematic-Based Approaches in Daxue Shan, Southeast Tibetan Plateau

: Active rock glaciers (ARGs) are important permafrost landforms in alpine regions. Identifying ARGs has mainly relied on visual interpretation of their geomorphic characteristics with optical remote sensing images, while mapping ARGs from their kinematic features has also become popular in recent years. However, a thorough comparison of geomorphic- and kinematic-based inventories of ARGs has not been carried out. In this study, we employed a multi-temporal interferometric synthetic aperture radar (InSAR) technique to derive the mean annual surface displacement velocity over the Daxue Shan, Southeast Tibet Plateau. We then compiled a rock glacier inventory by synergistically interpreting the InSAR-derived surface displacements and geomorphic features based on Google Earth images. Our InSAR-assist kinematic-based inventory (KBI) was further compared with a pre-existing geomorphic-based inventory (GBI) of rock glaciers in Daxue Shan. The results show that our InSAR-assist inventory consists of 344 ARGs, 36% (i.e., 125) more than that derived from the geomorphic-based method (i.e., 251). Only 32 ARGs in the GBI are not included in the KBI. Among the 219 ARGs detected by both approaches, the ones with area differences of more than 20% account for about 32% (i.e., 70 ARGs). The mean downslope velocities of ARGs calculated from InSAR are between 2.8 and 107.4 mm · a − 1 . Our comparative analyses show that ARGs mapping from the InSAR-based kinematic approach is more efﬁcient and accurate than the geomorphic-based approach. Nonetheless, the completeness of the InSAR-assist KBI is affected by the SAR data acquisition time, signal decorrelation, geometric distortion of SAR images, and the sensitivity of the InSAR measurement to ground deformation. We suggest that the kinematic-based approach should be utilized in future ARGs-based studies such as regional permafrost distribution assessment and water storage estimates.


Introduction
Rock glaciers are periglacial landforms developed in the cold mountain environment, consisting of ice and frozen rock debris that creep downslope due to the deformation of their internal ice-rock mixture [1,2]. They are visible indicators of permafrost existence in periglacial areas, and elevations at the front of rock glaciers are commonly referred to as the lower limit of discontinuous alpine permafrost [3]. Rock glaciers generally convey ice-rock mixture from the upper erosional part toward the lower depositional part (lobe) through an ill-defined transition zone (rooting zone) ( Figure 5) [4], which has important hydrological and geomorphological significance in alpine regions [5][6][7]. The creeping speeds of rock glaciers are sensitive to ground-thermal conditions. The response of rock glacier kinematics weaknesses of kinematic-and geomorphic-based approaches to compiling rock glaciers and aims at providing practical strategies for compiling rock glaciers completely and accurately. We selected the Daxue Shan, Southeast Tibetan Plateau, as the test area because a geomorphic-based rock glacier inventory in this region has been compiled [22]. Based on Sentinel-1A ascending SAR images acquired between 2015 and 2019, we derived a five-year-long LOS mean velocity map of the study area. We then visually examined the kinematic and geomorphic features of the active landforms and mapped the ARGs in the Daxue Shan region. We also compared our InSAR-derived inventory with the previous geomorphic-based rock glacier inventory [22]. We finally discuss the feasibilities of the geomorphic-and kinematic-based approach, respectively, and give suggestions on compiling a complete and accurate ARG inventory. We also discuss the possible influences of different approaches used to inventory rock glaciers on future rock glacier-based studies such as permafrost distribution assessment and water storage estimates.

Geological Overview
Daxue Shan is a mountain range in the Southeast Tibetan Plateau with an average elevation of about 4000 m a.s.l. It is also the geographical boundary of the Tibetan Plateau in the western Sichuan Basin, China ( Figure 1a). Our study area is located in the middle of Daxue Shan, between 29.95 • -30.58 • N and 101.51 • -101.95 • E, indicated by the red polygon in Figure 1a. This area is of great significance in regulating the regional climate, stabilizing the ecological environment in the middle and upper reaches of the Yangtze River, conserving water resources, and preventing soil erosion.
The climate in Daxue Shan is relatively warm and rich in rainfall because the region is one of the main passages of the Indian Ocean warm current entering the Chinese hinterland. The Daofu weather station (2957 m a.s.l) in the western slope of Daxue Shan documents 613.5 mm of the mean annual precipitation and a mean annual air temperature of 8.14 • C. The climatic condition in Daxue Shan is highly conducive to the formation of maritime glaciers and periglacial landforms [38]. Additionally, this region is located at the intersection of the Songpan, Chuandian, and South China tectonic blocks, which is a transition area with complex geological structures [39]. The major exposed strata in Daxue Shan [22] is Tertiary monzogranite, whose lithological and mineralogical characteristics are favorable for the development of rock glaciers [13]. The steep terrain, combined with intense tectonic activity, drives the sediment transport processes and the development of local rock glaciers [22].

Geomorphic-Based Rock Glacier Inventory in Daxue Shan
Ran and Liu [22] inventoried 295 rock glaciers in the Daxue Shan area through visual interpretation of geomorphic characteristics from Google Earth images. The geomorphic parameters of these rock glaciers were also documented with the support of ArcGIS software [22]. The GBI from Ran and Liu [22] showed that the rock glaciers in Daxue Shan cover a total area of about 55.7 km 2 , with altitudes concentrating between 4300 and 4600 m a.s.l. Based on the analysis of the distributions and environment factors of rock glaciers, Ran and Liu [22] suggested that the formation and development of rock glaciers are dominantly controlled by the regional topography, climate, and lithology.
Based on the genetic criteria, rock glaciers can be further divided into moraine type and talus type. The former type generally develops downstream from the end of existing glaciers or the cirque of degraded glaciers, while the latter type often develops below the steep rock wall that is covered with perennial snow [20,30]. Ran and Liu [22] reported that the moraine-type rock glaciers account for 49% (146) of their inventory, while the talus-type account for 51% (149).  [38] shown on the topographic map. The red polygon indicates the study area. (b) The distributions of ARGs mapped by kinematic-and geomorphic-based approaches. A circle filled with red or green represents the ARG mapped with a common boundary (Common Bound) or different boundaries (Diff. Bound) in the two inventories. The blue and yellow circles represent the ARGs only identified in the kinematic-and geomorphic-based inventories, respectively. The gray dots represent the relict rock glaciers inventoried from the geomorphic-based approach. (c) The outlines of ARGs from KBI shown on the mean annual surface velocity map derived from multi-temporal InSAR. The velocity map was divided into four subzones for a clear presentation, labeled as 1-4 and shown in Figure 4a-d.

Geomorphic-Based Rock Glacier Inventory in Daxue Shan
Ran and Liu [22] inventoried 295 rock glaciers in the Daxue Shan area through visual interpretation of geomorphic characteristics from Google Earth images. The geomorphic parameters of these rock glaciers were also documented with the support of ArcGIS software [22]. The GBI from Ran and Liu [22] showed that the rock glaciers in Daxue Shan cover a total area of about 55.7 km 2 , with altitudes concentrating between 4300 and 4600 m a.s.l. Based on the analysis of the distributions and environment factors of rock glaciers, Ran and Liu [22] suggested that the formation and development of rock glaciers are dominantly controlled by the regional topography, climate, and lithology.  [38] shown on the topographic map. The red polygon indicates the study area. (b) The distributions of ARGs mapped by kinematic-and geomorphic-based approaches. A circle filled with red or green represents the ARG mapped with a common boundary (Common Bound) or different boundaries (Diff. Bound) in the two inventories. The blue and yellow circles represent the ARGs only identified in the kinematicand geomorphic-based inventories, respectively. The gray dots represent the relict rock glaciers inventoried from the geomorphic-based approach. (c) The outlines of ARGs from KBI shown on the mean annual surface velocity map derived from multi-temporal InSAR. The velocity map was divided into four subzones for a clear presentation, labeled as 1-4 and shown in Figure 4a

Detecting Surface Deformation with Multi-Temporal InSAR
We employed the small-baseline subset (SBAS) multi-temporal InSAR method to derive the mean surface displacement velocity of the study area based on 116 ascending Sentinel-1A SAR images acquired between 2015 and 2019. The SBAS method works by exploring the SAR phase information on distributed scatterers, which is helpful for retaining Remote Sens. 2021, 13, 4931 5 of 22 SAR interferometric phase coherence and is suitable for mapping ground displacements in alpine regions [40].
The commercial GAMMA software was used to generate SAR interferograms with main steps of image coregistration, differential phase interferometry, phase filtering, unwrapping, and terrain geocoding [41]. We set the maximum temporal baseline and spatial baseline threshold of 84 days and 120 m, respectively. We also removed the interferograms with average spatial coherence smaller than 0.65 and finally obtained a fully connected network consisting of 532 interferograms ( Figure 2). The 30 m resolution AW3D30 DSM, provided by Japan Aerospace Exploration Agency (JAXA) [42], was used to circumvent the topographic phase. After generating SAR interferograms, the SAR phase was filtered using a Goldstein filter and unwrapped using a minimum cost flow algorithm. We then implemented the SBAS InSAR processing with open-source Mintpy software [43] to address mean surface displacement velocity. We used the EAR-5 products from the European Center for Medium-Range Weather Forecast to correct tropospheric phase delay [44].

Detecting Surface Deformation with Multi-Temporal InSAR
We employed the small-baseline subset (SBAS) multi-temporal InSAR method to derive the mean surface displacement velocity of the study area based on 116 ascending Sentinel-1A SAR images acquired between 2015 and 2019. The SBAS method works by exploring the SAR phase information on distributed scatterers, which is helpful for retaining SAR interferometric phase coherence and is suitable for mapping ground displacements in alpine regions [40].
The commercial GAMMA software was used to generate SAR interferograms with main steps of image coregistration, differential phase interferometry, phase filtering, unwrapping, and terrain geocoding [41]. We set the maximum temporal baseline and spatial baseline threshold of 84 days and 120 m, respectively. We also removed the interferograms with average spatial coherence smaller than 0.65 and finally obtained a fully connected network consisting of 532 interferograms ( Figure 2). The 30 m resolution AW3D30 DSM, provided by Japan Aerospace Exploration Agency (JAXA) [42], was used to circumvent the topographic phase. After generating SAR interferograms, the SAR phase was filtered using a Goldstein filter and unwrapped using a minimum cost flow algorithm. We then implemented the SBAS InSAR processing with open-source Mintpy software [43] to address mean surface displacement velocity. We used the EAR-5 products from the European Center for Medium-Range Weather Forecast to correct tropospheric phase delay [44].

Rock Glacier Mapping
With the mean ground displacement velocity map from SBAS InSAR processing available, we then mapped the ARGs in Daxue Shan using the workflow illustrated in Figure 3. First, we detected moving targets on the velocity map. A candidate moving target could be identified if at least 20-30 pixels homogenously deform [9]. We overlaid the velocity map on the optical satellite images and analyzed the geomorphic features of

Rock Glacier Mapping
With the mean ground displacement velocity map from SBAS InSAR processing available, we then mapped the ARGs in Daxue Shan using the workflow illustrated in Figure 3. First, we detected moving targets on the velocity map. A candidate moving target could be identified if at least 20-30 pixels homogenously deform [9]. We overlaid the velocity map on the optical satellite images and analyzed the geomorphic features of the moving targets to distinguish ARGs from other alpine landforms, such as glaciers, freezethaw mudflows, and landslides. ARGs generally develop in the lower reaches of alpine glaciers or lay at the flat section of valleys. Meanwhile, evidence of long-term downslope movements, such as transversal or longitudinal flow features, steep frontal and side slopes, and rare vegetation appearance, could be identified as distinct geomorphic features of ARGs. Here we used high-resolution Google Earth images for visual interpretation. Then, we delineated the outlines of ARGs from combined interpretations of their kinematic and geomorphic features. This was done by identifying the moving boundaries and geomorphic boundaries based on the kinematic and geomorphic features of the identified ARGs, respectively, followed by an integration of the two boundaries. Note that the interpreted boundary merely relying on the kinematic features may cause underestimation because it is common that only part of a rock glacier moves, while the interpretations only from the geomorphic features often lead to multiple results. We here chose the smallest geomorphic boundary greater than or equal to the moving boundary as the outline of an ARG.
boundaries and geomorphic boundaries based on the kinematic and geomorphic features of the identified ARGs, respectively, followed by an integration of the two boundaries. Note that the interpreted boundary merely relying on the kinematic features may cause underestimation because it is common that only part of a rock glacier moves, while the interpretations only from the geomorphic features often lead to multiple results. We here chose the smallest geomorphic boundary greater than or equal to the moving boundary as the outline of an ARG. We also calculated the geomorphic parameters (including ILP and FLP altitude, slope angle, aspect angle, area, and length) following the method proposed by Wang et al. [34]. First, the initial line point (ILP) and front line point (FLP), representing the places where permafrost creeping starts and ends of each ARG, were determined. For a detailed definition and identification method of these two points, see Section 3.3 of Wang et al. [34]. The 3D coordinates of ILPs and FLPs were extracted using the AW3D30 DSM, from We also calculated the geomorphic parameters (including ILP and FLP altitude, slope angle, aspect angle, area, and length) following the method proposed by Wang et al. [34]. First, the initial line point (ILP) and front line point (FLP), representing the places where permafrost creeping starts and ends of each ARG, were determined. For a detailed definition and identification method of these two points, see Section 3.3 of Wang et al. [34]. The 3D coordinates of ILPs and FLPs were extracted using the AW3D30 DSM, from which the lengths, slope angles (θ sl p ), and the aspect angles (α) of ARGs were further calculated.
In addition, we used Equation (1) to convert the velocities in the SAR LOS direction (V LOS ) into the downslope direction (V sl p ) [31]: where ϕ is the SAR flight azimuth angle and θ inc is the incidence angle. We calculated the mean velocity of all pixels within each ARG to represent its overall mobility. Considering the error propagation in Equation (1), we only calculated the velocities of ARGs whose conversion coefficients (i.e., the denominator part of Equation (1)) are higher than 0.2.

Comparison with the Existing Rock Glacier Inventory
To avoid biased comparison, we visually identified the ARGs from GBI made by Ran and Liu (2018) and also calculated the geomorphic parameters using the method mentioned above. We then compared the distributions, outlines, and geomorphic parameters of ARGs from the GBI with our InSAR-assist kinematic-based inventory.
Note that the GBI from Ran and Liu [22] did not contain the activity status of the rock glaciers. To compare GBI with our inventory, we identified the active and transitional rock glaciers in the GBI following the guidelines formulated by the IPA Action Group [9]. Among the 295 rock glaciers from GBI, 251 rock glaciers were classified into ARGs (including both the active and transitional ones), and the others (44) are relict rock glaciers (Figure 1b). Results were summarized into four cases: ARGs sharing the common boundary in KBI and GBI, ARGs having different outlines in KBI and GBI, ARGs only identified in KBI, and ARGs only identified in GBI (Figure 1b). Considering the relative subjectivity in mapping rock glacier outlines, the so-called common boundary of an ARG means that the polygons from two inventories have the same geomorphological structure, and the area difference is less than 20%.

ARGs from the InSAR-Assist Kinematic-Based Inventorying
We mapped a total of 344 ARGs through the InSAR-assist kinematic-based approach in the Daxue Shan area. Figure 4 shows the outlines of ARGs from KBI overlaying on the ground surface velocity map derived from SBAS InSAR. The ARGs in the Daxue Shan region are widely distributed in the elevation range 3860-5094 m a.s.l., among which the moraine-type and talus-type ARGs account for 52% (180) and 48% (164), respectively. The ARGs mainly face the northeast.

Comparisons of ARG Distributions and Outlines in KBI and GBI
Two hundred nineteen ARGs were included in both the GBI of Ran and Liu [22] and our InSAR-assist KBI. We found that 149 ARGs from KBI and GBI (accounting for 43% of KBI and 60% of GBI, respectively) were mapped with common boundaries. Figure 5a-d shows the examples of two ARGs sharing common boundaries. The ARG1 in Figure 5a-b is a composite ARG, and the sediment is supplied from two hillslopes. The ARG2 in Figure 5c-d is a typical tongue-shaped ARG, consisting of one rooting zone and one lobe. All the rooting zones and lobes of these two ARGs were mapped in both of the two inventories, with area differences less than 20%. The InSAR-assist KBI newly identified 125 ARGs (accounting for 36% of the total number of KBI). Figure 6a,b shows two examples of the newly identified ARGs. The geomorphic features indicating permafrost-related creeping of ARG5 and ARG6 are unclear, making them difficult to distinguish from the surrounding ground in the optical image ( Figure 6a). In contrast, evident moving signals of these two ARGs are presented in the multi-temporal InSAR displacement velocity map ( Figure 6b).
Thirty-two ARGs were included in GBI but not identified in KBI (accounting for 13% of GBI). We found that about half of these ARGs (15) have aspects of the near north-south direction, indicating that these ARGs are omitted probably because the InSAR measurement is not sensitive to the ground movements in this direction. As an example, Figure 6c-d shows an ARG7 mainly flowing toward the SAR flight azimuth direction (i.e., near north-south direction), while little deformation signal is observed in the InSAR displacement velocity map. We also found 17 ARGs that do not creep along the north-south direction but show geomorphic evidence of permafrost-related creeping (see ARG8 Figure 6e-f as an example). This indicates that the InSAR-assist KBI also has its limitations, which are discussed in Section 5.1.  Table 1 shows the overall statistics of the GBI and KBI. The largest ARGs in GBI and KBI have areas of 1.012 and 0.92 km 2 , respectively. The ARGs in KBI are generally shorter and steeper in slope than the GBI. The ARGs in KBI and GBI have similar altitude distributions, with the former showing an altitude range of 3860-5094 m a.s.l., and the latter showing a distribution range of 3823-5095 m a.s.l. The total area of the ARGs in KBI is 24% larger than that of the GBI. Next, we show in detail the differences in the geomorphic parameters between GBI and KBI.

Comparisons of ARG Distributions and Outlines in KBI and GBI
Two hundred nineteen ARGs were included in both the GBI of Ran and Liu [22] and our InSAR-assist KBI. We found that 149 ARGs from KBI and GBI (accounting for 43% of KBI and 60% of GBI, respectively) were mapped with common boundaries. shows the examples of two ARGs sharing common boundaries. The ARG1 in Figure 5a,b is a composite ARG, and the sediment is supplied from two hillslopes. The ARG2 in Figure 5c,d is a typical tongue-shaped ARG, consisting of one rooting zone and one lobe. All the rooting zones and lobes of these two ARGs were mapped in both of the two inventories, with area differences less than 20%.
We found 70 ARGs (accounting for 20% of KBI and 28% of GBI) having different outlines. Figure 5e,f shows two adjacent ARGs as examples that have different outlines in the two inventories. The sediments at both ARGs are supplied from multiple headwalls and transferred into multiple lobes. Two rooting zones and two lobes of the ARG3 were excluded in the GBI. All the rooting zones of ARG4 were partially mapped, with one lobe excluded in the GBI. The areas of ARG3 and ARG4 mapped in KBI are 51% and 29% larger than those from the GBI. It can be seen that, for these complex situations as illustrated above, accurate delineation of the ARG outlines is challenging by only interpreting the geomorphic features. The uncertainties of outlining can be reduced to a large extent by using the InSAR-assist kinematic approach.
The InSAR-assist KBI newly identified 125 ARGs (accounting for 36% of the total number of KBI). Figure 6a,b shows two examples of the newly identified ARGs. The geomorphic features indicating permafrost-related creeping of ARG5 and ARG6 are unclear, making them difficult to distinguish from the surrounding ground in the optical image ( Figure 6a). In contrast, evident moving signals of these two ARGs are presented in the multi-temporal InSAR displacement velocity map (Figure 6b).
Thirty-two ARGs were included in GBI but not identified in KBI (accounting for 13% of GBI). We found that about half of these ARGs (15) have aspects of the near north-south direction, indicating that these ARGs are omitted probably because the InSAR measurement is not sensitive to the ground movements in this direction. As an example, Figure 6c,d shows an ARG7 mainly flowing toward the SAR flight azimuth direction (i.e., near northsouth direction), while little deformation signal is observed in the InSAR displacement velocity map. We also found 17 ARGs that do not creep along the north-south direction but show geomorphic evidence of permafrost-related creeping (see ARG8 Figure 6e,f as an example). This indicates that the InSAR-assist KBI also has its limitations, which are discussed in Section 5.1. Table 1 shows the overall statistics of the GBI and KBI. The largest ARGs in GBI and KBI have areas of 1.012 and 0.92 km 2 , respectively. The ARGs in KBI are generally shorter and steeper in slope than the GBI. The ARGs in KBI and GBI have similar altitude distributions, with the former showing an altitude range of 3860-5094 m a.s.l., and the latter showing a distribution range of 3823-5095 m a.s.l. The total area of the ARGs in KBI is 24% larger than that of the GBI. Next, we show in detail the differences in the geomorphic parameters between GBI and KBI. Table 1. Statistical summary of ARGs inventoried in KBI and GBI. Note K-Moraine, K-Talus, and K-All represent the moraine-type, talus-type, and all ARGs from KBI. Similarly, G-Moraine, G-Talus, and G-All represent the moraine-type, talus-type, and all ARGs from GBI. The total area of ARGs in KBI is 66.6 km 2 , among which the moraine type accounts for 58% (38.7 km 2 ), and the talus type accounts for 42% (27.9 km 2 ). In contrast, the total area of ARGs in GBI is 50.6 km 2 , among which the moraine type accounts for 53% (26.9 km 2 ), and the talus type accounts for 47% (23.8 km 2 ) ( Table 1). The mean area of all the ARGs from GBI is 0.2, slightly higher than that from KBI (0.19). Most of the ARGs in both KBI and GBI have areas less than 0.3 km 2 (Figure 7a,b).  The mean length of the moraine-type, talus-type, and all ARGs in GBI are 839, 840, and 839 m, overall longer than those of KBI (792, 682, and 736 m, respectively) ( Table 1); 62% and 60% of the ARGs in KBI and GBI have lengths ranging from 400 to 1000 m (Figure 7c,d).

Comparative Analysis of the Geomorphic Parameters in GBI and KBI
Figure 7e,f shows that the ARGs from GBI generally have gentler slopes than those from the KBI. Additionally, the slopes of ARGs from KBI are more widely distributed than GBI. The slopes of ARGs in KBI mainly range between 13 • and 25 • , while the slopes in GBI mainly range between 11 • and 20 • (Figure 7e). Particularly, the ARGs newly identified in KBI mostly have slopes higher than 20 • , probably because the steep terrain hinders correct interpretations of the geomorphic features and leads to the omissions of ARGs in GBI. cause moraine-type ARGs are generally developed from ablated glaciers or interact with glaciers and ice patches within a confined altitude range. We can also observe that the talus-type ARGs reach lower altitudes than moraine-type ARGs in KBI, while this pattern is not obvious in GBI (Figure 8c).

Aspect Angle
Both the GBI and KBI have the most significant number of ARGs (respectively accounting for 24% and 28%) facing the northeast, indicating the northeast slopes are favorable for the ARG formation in Daxue Shan. Except for the aggregated distribution in the northeast direction of GBI, a small amount of ARGs in the GBI face the south (4%) and southeast (4%), and the number of ARGs facing the other directions are nearly equal (11%-15% in each direction). In KBI, ARGs show aggregated distribution in the northeast direction (28%), east direction (21%), and north direction (15%), and are also almost evenly distributed in the other directions (7%-9% in each direction) (Figure 9a). The significant NE-weighted distribution of ARGs likely results from the topographic control of the Daxue Shan. Previous studies have shown that the slopes receiving lower solar radiation generally favor the presence of ARGs because this helps to preserve ground ice inside ARGs [34]. Given the strike angle of near NW15° for the Daxue Shan mountain ridge, the slopes facing northeast are thus the favorable slopes for forming ARGs.
The moraine-type ARGs in both KBI and GBI show similar aspect distribution pat- The distribution ranges of ILP altitudes of talus-type ARGs are obviously larger than that of the moraine-type ARGs in both inventories, indicating that moraine-type ARGs are generally distributed in a narrow altitudinal range (Figure 8a,c). This is probably because moraine-type ARGs are generally developed from ablated glaciers or interact with glaciers and ice patches within a confined altitude range. We can also observe that the talus-type ARGs reach lower altitudes than moraine-type ARGs in KBI, while this pattern is not obvious in GBI (Figure 8c).

Aspect Angle
Both the GBI and KBI have the most significant number of ARGs (respectively accounting for 24% and 28%) facing the northeast, indicating the northeast slopes are favorable for the ARG formation in Daxue Shan. Except for the aggregated distribution in the northeast direction of GBI, a small amount of ARGs in the GBI face the south (4%) and southeast (4%), and the number of ARGs facing the other directions are nearly equal (11-15% in each direction). In KBI, ARGs show aggregated distribution in the northeast direction (28%), east direction (21%), and north direction (15%), and are also almost evenly distributed in the other directions (7-9% in each direction) (Figure 9a). The significant NE-weighted distribution of ARGs likely results from the topographic control of the Daxue Shan. Previous studies have shown that the slopes receiving lower solar radiation generally favor the presence of ARGs because this helps to preserve ground ice inside ARGs [34]. Given the strike angle of near NW15 • for the Daxue Shan mountain ridge, the slopes facing northeast are thus the favorable slopes for forming ARGs.

Active Rock Glacier Velocities
Compared to GBI, the most prominent advantage of KBI is that it can quantify the surface motions of ARGs directly. We calculated the mean downslope movement velocities of 312 ARGs (including 167 moraine-type ARGs and 145 talus-type ARGs) following Equation (1). The velocities of the other 32 ARGs were not documented because their low LOS-to-downslope conversion coefficient may introduce large uncertainties (see Section 3.2). The downslope movement velocities of ARGs in Daxue Shan show high spatial heterogeneity. The mean annual velocities of ARGs range from 2.8 to 107.4 mm•a -1 , with a mean velocity of about 25 mm•a -1 for all the ARGs ( Table 2). The mean movement velocity of moraine-type ARGs (25 mm•a -1 ) is slightly higher than that of talus-type ARGs (23.7 mm•a -1 ) ( Table 2).  Figure 10 shows the violin plot and histogram distribution of the downslope movement velocities for the mapped ARGs. Most velocities concentrate between 10 and The moraine-type ARGs in both KBI and GBI show similar aspect distribution patterns (Figure 9b), while the talus-type ARGs in the two inventories show different distributions. Seventy-six percent of the talus-type ARGs in KBI face the north, northeast, and east, whereas only 37% of the talus-type ARGs in GBI face these directions.

Active Rock Glacier Velocities
Compared to GBI, the most prominent advantage of KBI is that it can quantify the surface motions of ARGs directly. We calculated the mean downslope movement velocities of 312 ARGs (including 167 moraine-type ARGs and 145 talus-type ARGs) following Equation (1). The velocities of the other 32 ARGs were not documented because their low LOS-to-downslope conversion coefficient may introduce large uncertainties (see Section 3.2). The downslope movement velocities of ARGs in Daxue Shan show high spatial heterogeneity. The mean annual velocities of ARGs range from 2.8 to 107.4 mm·a −1 , with a mean velocity of about 25 mm·a −1 for all the ARGs ( Table 2). The mean movement velocity of moraine-type ARGs (25 mm·a −1 ) is slightly higher than that of talus-type ARGs (23.7 mm·a −1 ) ( Table 2).  Figure 10 shows the violin plot and histogram distribution of the downslope movement velocities for the mapped ARGs. Most velocities concentrate between 10 and 35 mm·a −1 (Figure 10a). The velocities for both moraine-type and talus-type ARGs mainly range between 10 and 50 mm·a −1 . Only two ARGs have creeping speeds higher than 100 mm·a −1 , and both are moraine-type. The velocity range of ARGs in Daxue Shan is comparable to that in the Uinta Mountains reported by Brencher, Handwerger and Munroe [36]. However, the velocity range is lower than other alpine regions in the Alps, Sierra Nevada, and Northern Tienshan [26,31,34,35].

The Geomorphic-Based and Kinematic-Based ARG Inventory Approaches
The geomorphic-based approach can compile exhaustive rock glaciers, including moving and non-moving ones, by visual interpretation of aerial images [9], without tedious data processing. Rock glaciers from GBI have been frequently used to infer the linkages between rock glacier distributions and regional environment variables, e.g., [13,22,45]. However, the completeness of rock glacier mapping and the correctness of its activity classification through the geomorphic-based approach highly depend on the image quality and expertise of the operssator, hence always leading to somewhat subjective mapping results [4]. For example, 36% of ARGs identified by the InSAR-assist kinematic-based approach were failed to be included in GBI in Daxue Shan. Additionally, 29% of ARGs from GBI were mapped incompatibly with respect to their kinematic boundaries. These variabilities in mapping outcomes can propagate to the subsequent geomorphic parameter calculations, and even induce large uncertainties in the further cryosphere-related studies relying on rock glaciers Compared to the geomorphic-based approach, the kinematic-based approach can help effectively detect and accurately delineate ARGs. The kinematic-based approach can also quantify the surface motions of ARGs. This approach is thus relatively more objective and can reduce the uncertainties in future studies based on the identified ARGs [46]. We here suggest using multi-annual mean displacement signals to detect moving targets firstly, and then systematically integrating the kinematic features and geomorphic characteristics to identify and map rock glaciers.
Our study in Daxue Shan demonstrates that the InSAR-assist kinematic-based approach also has its own limitations. First, the ARGs facing the satellite flight direction may be missed by this approach because of the insensitivity of InSAR measurement to ground movement along this direction. In our study region, 15 ARGs facing the near north-south direction (accounting for 6% of GBI) failed to be detected by InSAR, alt-

The Geomorphic-Based and Kinematic-Based ARG Inventory Approaches
The geomorphic-based approach can compile exhaustive rock glaciers, including moving and non-moving ones, by visual interpretation of aerial images [9], without tedious data processing. Rock glaciers from GBI have been frequently used to infer the linkages between rock glacier distributions and regional environment variables, e.g., [13,22,45]. However, the completeness of rock glacier mapping and the correctness of its activity classification through the geomorphic-based approach highly depend on the image quality and expertise of the operssator, hence always leading to somewhat subjective mapping results [4]. For example, 36% of ARGs identified by the InSAR-assist kinematic-based approach were failed to be included in GBI in Daxue Shan. Additionally, 29% of ARGs from GBI were mapped incompatibly with respect to their kinematic boundaries. These variabilities in mapping outcomes can propagate to the subsequent geomorphic parameter calculations, and even induce large uncertainties in the further cryosphere-related studies relying on rock glaciers Compared to the geomorphic-based approach, the kinematic-based approach can help effectively detect and accurately delineate ARGs. The kinematic-based approach can also quantify the surface motions of ARGs. This approach is thus relatively more objective and can reduce the uncertainties in future studies based on the identified ARGs [46]. We here suggest using multi-annual mean displacement signals to detect moving targets firstly, and then systematically integrating the kinematic features and geomorphic characteristics to identify and map rock glaciers.
Our study in Daxue Shan demonstrates that the InSAR-assist kinematic-based approach also has its own limitations. First, the ARGs facing the satellite flight direction may be missed by this approach because of the insensitivity of InSAR measurement to ground movement along this direction. In our study region, 15 ARGs facing the near north-south direction (accounting for 6% of GBI) failed to be detected by InSAR, although all of them exhibit clear geomorphic features indicating permafrost-related creeping. Second, shadow and layover distortions in SAR images can degrade the displacement signal and lead to some of the ARGs being omitted in the inventory. Some researchers only included ARGs located outside layover and shadow zones in the ARG inventories [31]. The SAR images used in this study have about 5.9% pixels affected by overlay and shadow effects, implying that these SAR geometric distortions have limited influences on the inventory results. However, for regions with rugged topography where SAR images may be strongly affected by shadow and layover, it is suggested to use SAR images acquired from different orbits (e.g., ascending and descending) to conduct the inventory. In addition, SAR signal decorrelation due to the alpine environment changes and atmospheric delays is also the common error source in InSAR displacement measurements, which may influence the completeness and correctness of the kinematic-based rock glacier inventory.
Future ARG inventories should value the integration of both geomorphic and kinematic features of ARGs using advanced methods such as machine learning or deep learning algorithms. These advanced methods can automatically recognize recurring patterns and textures and thus reduce the workload in compiling rock glaciers. Furthermore, recent studies have shown that machine learning-based SAR data processing can better identify and correct the artifacts in surface displacement measurement [47]. We suggest that developing machine learning-based methods for extracting both geomorphic and kinematic information of rock glaciers is a possible solution to homogeneous compilation and global assemblage of rock glacier inventories.

Influences of SAR Image Selection on ARG Inventory
We discussed the advantages and limitations of InSAR-assist ARG mapping in Section 5.1. Besides using multi-temporal SAR images, as shown in this study, using a signal SAR image pair has also been frequently applied to map ARGs. Additionally, considering that the downslope movements of ARGs generally show a seasonal fluctuation pattern, we postulate that using SAR images acquired in different seasons may also influence the inventory. Here, we discuss the reasonable strategy of SAR data processing aimed at overcoming the limitations of InSAR for ARG inventory.

Active Rock Glacier Inventory from a Single SAR Image Pair
The use of only one or a few SAR interferograms to detect ARGs has been successfully implemented in the Sierra Nevada [31], Tien Shan [34], and Dry Andes Mountains [48], etc., where 59, 261, and 2116 ARGs have been mapped, respectively. However, the effectiveness of this strategy is closely related to the quality of the SAR interferogram. Given that ARGs creep with diverse velocities, SAR interferograms should be carefully selected with proper temporal baselines. For instance, ground displacements in the range 8-90 cm a −1 can be measured from a Sentinel-1 SAR interferogram with 12 day temporal baseline, while the displacement range will be reduced to 2-30 cm a −1 when the temporal baseline extends to one month [49]. The simple calculations imply that slowing ARGs may be missed by using a SAR interferogram with a temporal baseline that is too short because the ARGs cannot be distinguished from the stable surrounding ground. On the other hand, if the temporal baseline is too long, SAR temporal decorrelation may, in turn, cause unsuccessful measurement of the displacement signals. Selecting the SAR image pair with an appropriate temporal baseline is thus crucial for single SAR interferogram-based ARG inventorying.
To demonstrate the influence of interferogram time span on the variability of ARG mapping, we selected a region (see the location in Figure 4b) containing 11 ARGs and then examined the ARG numbers identified from three SAR interferograms with 12-, 36-, and 48-day temporal baselines (Figure 11c-e). We can identify 6, 10, and 4 ARGs (the black polygons) from the 12-, 36-, and 38-day temporal-span interferograms, respectively. For example, the ARG No. 3 (Figure 11a) cannot be identified in the 12-day-span interferogram, but show some movement signals in the 36-and 48-day-span interferogram. Additionally, ARG No. 10 cannot be identified in all three interferograms, but can be clearly detected in the multi-temporal InSAR mean velocity map (Figure 11b). This indicates that using a single SAR interferogram for ARG inventory could lead to a high omission rate.
Remote Sens. 2021, 13, x FOR PEER REVIEW 17 of 23 urement accuracy [26,31,34]. Nevertheless, our study in Daxue Shan suggests that using multi-temporal InSAR is a better option than using only one single SAR image pair for ARGs inventorying. Surface displacement velocity maps derived from the advanced InSAR algorithms can provide better kinematic information of ARGs and ensure completeness and accuracy of the inventory. Figure 11. Active rock glaciers in zone "Z" shown in Figure 4b, overlaid on the (a) Google Earth image and (b) multi-temporal InSAR velocity map; (c-e) three interferograms with 12-, 36-, and 48-day time span, respectively. The ARGs with black outlines represent that they can be identified effectively, while the ones with white outlines represent that they cannot be effectively detected.

Influence of SAR Acquisition Season on ARG inventories
The winter season in Daxue Shan is defined from October to March of the following year, and the summer season is defined from April to September within a year. To infer the influence of SAR acquisition season on ARG inventories, we derived the surface displacement velocities using images only acquired either in the winter (Figure 12a) or summer (Figure 12b) season, respectively. As a comparison, Figure 12c shows the velocity map using all the SAR images.
Both the winter and summer surface velocities are different from the full-data-derived velocities, especially in the mountain ridges and valleys (see the polygons in Figure 12). We randomly selected an area P as an example and compared the displacements for this region in the three velocity maps. The velocity values are positive for area P in the winter velocity map representing displacements toward the satellite in the LOS direction (Figure 12a). However, the velocity values become negative in the summer velocity map (Figure 12b). Widespread homogenous "displacement signals" can be found in the area P either in the winter or summer seasons, preventing us from iden- Figure 11. Active rock glaciers in zone "Z" shown in Figure 4b, overlaid on the (a) Google Earth image and (b) multi-temporal InSAR velocity map; (c-e) three interferograms with 12-, 36-, and 48-day time span, respectively. The ARGs with black outlines represent that they can be identified effectively, while the ones with white outlines represent that they cannot be effectively detected.
When only using one image pair to inventory ARGs, to improve the compilation completeness, it is suggested to use SAR images acquired by the SAR sensor emitting long radar wavelengths, such as the L-band ALOS PALSAR. Long radar wavelength generally has a good penetration ability and thus can preserve SAR coherence for measurement accuracy [26,31,34]. Nevertheless, our study in Daxue Shan suggests that using multitemporal InSAR is a better option than using only one single SAR image pair for ARGs inventorying. Surface displacement velocity maps derived from the advanced InSAR algorithms can provide better kinematic information of ARGs and ensure completeness and accuracy of the inventory.

Influence of SAR Acquisition Season on ARG Inventories
The winter season in Daxue Shan is defined from October to March of the following year, and the summer season is defined from April to September within a year. To infer the influence of SAR acquisition season on ARG inventories, we derived the surface displacement velocities using images only acquired either in the winter (Figure 12a) or summer (Figure 12b) season, respectively. As a comparison, Figure 12c shows the velocity map using all the SAR images.
Remote Sens. 2021, 13, x FOR PEER REVIEW 18 of 23 tifying the moving ground targets. This indicates that using the ground deformation products generated from only winter or summer SAR images can result in a highly biased inventory. The inverse displacement signals for the same area between the winter and summer velocity map probably result from seasonal atmospheric delay or hydrological loading effects [50]. In contrast, clear displacement signals associated with ARG creeping can be seen in the velocity map derived from full SAR data (Figure 12c), suggesting that the SAR dataset used for ARG inventories should cover at least one complete cycle year.

Implication for Rock Glacier-Based Permafrost Distribution Estimation
Rock glaciers have been ubiquitously used in permafrost distribution mapping since they are regarded as the visible indicators of permafrost occurrence. Based on rock glaciers' activity status (i.e., active or relict), the probability of permafrost occurrence can be linked to the regional environmental variables such as altitude, mean annual air temperature (MMAT), ground temperature, and solar radiation through statistic models. A few studies have relied on the geomorphic-based rock glacier inventories to estimate the permafrost distribution in the Southern Alps, New Zealand [51], Hunza River Basin, Chilean Andes [52], Western Karakoram [53], and Sikkim Himalayas [54].
We compared the ARGs mapped from GBI and KBI with the existing permafrost distribution maps from statistical modeling. Figure 13a,b shows the GBI and KBI superimposed on the Global Permafrost Zonation Index (PZI) map developed by Gruber [55] and Northern Hemisphere Permafrost Occurrence Probability (POP) map developed by Obu, et al. [56]. The PZI map was established by modeling the relationship between long-term mean annual air temperature (MAAT) and permafrost occurrence using the Both the winter and summer surface velocities are different from the full-data-derived velocities, especially in the mountain ridges and valleys (see the polygons in Figure 12). We randomly selected an area P as an example and compared the displacements for this region in the three velocity maps. The velocity values are positive for area P in the winter velocity map representing displacements toward the satellite in the LOS direction (Figure 12a). However, the velocity values become negative in the summer velocity map (Figure 12b). Widespread homogenous "displacement signals" can be found in the area P either in the winter or summer seasons, preventing us from identifying the moving ground targets. This indicates that using the ground deformation products generated from only winter or summer SAR images can result in a highly biased inventory. The inverse displacement signals for the same area between the winter and summer velocity map probably result from seasonal atmospheric delay or hydrological loading effects [50]. In contrast, clear displacement signals associated with ARG creeping can be seen in the velocity map derived from full SAR data (Figure 12c), suggesting that the SAR dataset used for ARG inventories should cover at least one complete cycle year.

Implication for Rock Glacier-Based Permafrost Distribution Estimation
Rock glaciers have been ubiquitously used in permafrost distribution mapping since they are regarded as the visible indicators of permafrost occurrence. Based on rock glaciers' activity status (i.e., active or relict), the probability of permafrost occurrence can be linked to the regional environmental variables such as altitude, mean annual air temperature (MMAT), ground temperature, and solar radiation through statistic models. A few studies have relied on the geomorphic-based rock glacier inventories to estimate the permafrost distribution in the Southern Alps, New Zealand [51], Hunza River Basin, Chilean Andes [52], Western Karakoram [53], and Sikkim Himalayas [54].
We compared the ARGs mapped from GBI and KBI with the existing permafrost distribution maps from statistical modeling. Figure 13a,b shows the GBI and KBI superimposed on the Global Permafrost Zonation Index (PZI) map developed by Gruber [55] and Northern Hemisphere Permafrost Occurrence Probability (POP) map developed by Obu, et al. [56]. The PZI map was established by modeling the relationship between longterm mean annual air temperature (MAAT) and permafrost occurrence using the simple but widely accepted cumulative normal distribution [55]. The POP map used the TTOP (temperature at the top of permafrost) equilibrium model to link permafrost occurrence with environment variables, including remotely sensed land surface temperatures and downscaled ERA-Interim climate reanalysis data, tundra wetness classes, and landcover map [57]. Both the PZI and POP values close to 1 indicate high probabilities of permafrost occurrence. Statistics show that the mean PZI values are both smaller than 0.1 for the ARGs from KBI (0.07) and GBI (0.09) in Daxue Shan. The mean POP values for ARGs in KBI and GBI are moderate, with values of 0.48 and 0.46, respectively (Table 3). Among the 344 and 251 ARGs included in KBI and GBI, only 2 and 1 ARGs have PZI values higher than 0.5. The numbers of ARGs showing POP values higher than 0.5 are 164 and 111 for KBI and GBI, respectively.
The KBI and GBI have comparable permafrost statistical indices with slight discrepancy (Table 3), although our results (Section 4) have shown that the ARGs from these two inventories are different in distributions, outlines, and geographic parameters. The similar statistics between PZI and POP of GBI and KBI are probably due to the coarse spatial resolution (both at 1 km 2 ) of PZI and POP maps.
The small to moderate PZI and POP values on the ARGs also indicate that the global permafrost distribution products from statistical modeling significantly underestimated the permafrost extent in the Daxue Shan region. This is probably due to the lack of accurate constraints (i.e., air temperature and precipitation) on the permafrost estimation model. Furthermore, given that solar radiation, topography, and snow depth all have influences on the ground thermal regime [58], which control the occurrence of permafrost, using relatively simple empirical models to estimate large-scale permafrost distribution may somewhat be biased. In addition, ARGs move downslope and may reach the permafrostfree region. Due to the cooling effect of rock glaciers through their high porosity and storage of ice, rock glaciers can lead to a more favorable environment for permafrost occurrence than the surrounding terrain [58,59].   The small to moderate PZI and POP values on the ARGs also indicate that the global permafrost distribution products from statistical modeling significantly underestimated

Implications for Water Storage Content Estimation
Although ARGs are generally smaller than glaciers, they are regarded as long-term reservoirs of great value in hydrology because of their large number and wide distribution. With increasing temperature in the future, ARGs may become increasingly more important freshwater resources in arid and semiarid areas [5]. Rock glaciers compiled from the geomorphic-based approach have been used for water content estimations in several mountain ranges, with assumptions of a certain ice content and ice thickness. For example, based on the 669 rock glaciers from geomorphic-based interpretation, Janke, Ng and Bellisario [17] assumed a rock glacier ice content of about 25-44.9% and estimated water storage of 2200-4858 km 3 in the semiarid Aconcagua River Basin of Chile. Similarly, Jones, et al. [60] estimated the water content of rock glaciers at about 1.672 × 10 10 -2.508 × 10 10 km 3 in Nepalese Himalaya with an assumption of ice content at 40-60% from a GBI including 1137 ARGs. Considering the heterogeneous properties of rock glaciers, apparently, the simple assumptions of ice content and thickness may lead to significantly biased estimates.
Using rock glaciers from KBI may overcome these limitations of GBI for regional water storage content calculations. First, as shown in the Daxue Shan case, the KBI can achieve a more complete compilation of rock glaciers than GBI, with 36% differences in numbers and 24% differences in total areas. Second, the kinematic information can help improve the accuracy of ice content estimates. Several studies have shown that rock glacier downslope movements can be related to their ice content through rheological models [12,61]. The better estimation of ice content can thus improve the water storage content estimates. We suggest that the kinematic-based rock glacier inventories are useful for inferring the hydrological significance of rock glaciers [61].

Conclusions
In this study, we compared the kinematic-based and geomorphic-based rock glacier inventory approaches in the Daxue Shan region. We employed the multi-temporal InSAR velocity map and Google Earth optical images to assist the identification of ARGs. Compared to the inventory from the geomorphic-based approach, our InSAR-assist kinematicbased inventory overall shows better completeness and accuracy. We newly identified 125 ARGs in the study region, which is 36% more than the GBI.
We found that 149 ARGs share common outlines and 70 ARGs were mapped with different polygon outlines. Thirty-two ARGs included in GBI were not identified by the InSAR-assist kinematic-based approach (accounting for 11% of GBI), which was mainly caused by the insensitivity of InSAR measurements to ground movements parallel to the satellite motion. The total area of ARGs in KBI is 24% more than that in GBI. The ARGs in GBI are overall longer and distributed on more gentle slopes than those from KBI. Both ILP and FLP altitudes in the two inventories show similar altitudinal distributions. Aspect distribution patterns of GBI and KBI show obvious discrepancies, especially for the talustype ones. The downslope velocities of ARGs are concentrated in the range 10-35 mm·a −1 , with only 2 ARGs creeping faster than 100 mm·a −1 .
Inventories from InSAR-assist kinematic-based approach may be affected by the season of SAR data acquisition, geometry distortions, the sensitivity of InSAR measurement to ground deformation, and SAR signal decorrelations. Using multi-temporal SAR images spanning at least one year to derive the surface deformation map is an effective way to overcome these limitations. We inferred that the inventory approach greatly influences the additional rock glacier-based studies, such as permafrost distribution mapping and water equivalent estimates. The kinematic-based inventory approach is useful in the future because of its potential in providing complete and accurate ARG inventories and realistic constraints on the inversion of hydrological parameters on ARGs. Advanced methods such as machine learning or deep learning algorithms for automatic identification of both geomorphic and kinematic features of ARGs should also be valued in future inventories.