Inventory and changes of rock glacier creep speeds in Ile Alatau and Kungöy Ala-Too, northern Tien Shan, since the 1950s

Spatio-temporal patterns related to the viscous creep in perennially frozen sediments of rock glaciers in cold mountains have rarely been studied outside the densely populated European Alps. This study investigates the spatial and temporal variability of rock glacier movement in the Ile Alatau and Kungöy Ala-Too mountain ranges, northern Tien Shan, a region with particularly large and fast rock glaciers. Over the study region of more than 3000 km, an inventory of slope movements is constructed using a large number of radar interferograms and high-resolution optical imagery. The inventory 20 includes more than 900 landforms, of which around 550 are interpreted as rock glaciers. Out of the active rock glaciers inventoried, 45 are characterised by a rate of motion exceeding 100 cm/a. From these fast rock glaciers we select six (Gorodetzky, Morenny, Archaly, Ordzhonikidze, Karakoram and Kugalan Tash) and study them in more detail using offset tracking between repeat airphotos, and historical and modern high resolution optical satellite data. Two of these rock glaciers show a steady increase of decadal surface velocities from the 1950s onwards with speeds being roughly two to four times 25 higher in recent years compared to the 1950s and 1960s. Three rock glaciers show similar accelerations over the recent one to two decades, but also phases of increased speeds in the 1960s. This development indicates a possible significant increase in current sediment and ice fluxes through rock glaciers and implies that their material transport in the region might gain geomorphodynamic importance relative to material transport by glaciers, assuming the latter decreases together with the regional glacier shrinkage. The study demonstrates how air and satellite image archives are exploited to construct one of the 30 longest decennial times series of rock glacier speeds currently available. Our results are in line with findings from Europe about rock glacier speeds increasing with atmospheric warming, and underline local variability of such overall response.

inventoried, 45 are characterised by a rate of motion exceeding 100 cm/a. From these fast rock glaciers we select six (Gorodetzky, Morenny, Archaly, Ordzhonikidze, Karakoram and Kugalan Tash) and study them in more detail using offset tracking between repeat airphotos, and historical and modern high resolution optical satellite data. Two of these rock glaciers show a steady increase of decadal surface velocities from the 1950s onwards with speeds being roughly two to four times 25 higher in recent years compared to the 1950s and 1960s. Three rock glaciers show similar accelerations over the recent one to two decades, but also phases of increased speeds in the 1960s. This development indicates a possible significant increase in current sediment and ice fluxes through rock glaciers and implies that their material transport in the region might gain geomorphodynamic importance relative to material transport by glaciers, assuming the latter decreases together with the regional glacier shrinkage. The study demonstrates how air and satellite image archives are exploited to construct one of the 30 longest decennial times series of rock glacier speeds currently available. Our results are in line with findings from Europe about rock glacier speeds increasing with atmospheric warming, and underline local variability of such overall response.

Introduction
Debris-ice mixtures on mountain slopes can accumulate, deform under gravity, and form striking tongue-shaped features of viscous flow, up to a kilometre wide and several kilometres longso called rock glaciers (or rockglaciers) (Haeberli, 1985;Martin and Whalley, 1987;Barsch, 1996). Rock glaciers can be defined as "the visible expression of cumulative deformation by long-term creep of ice/debris mixtures under permafrost conditions" (Berthling, 2011), a definition that we will follow in 5 this study. Prominent rock glaciers are found in most cold mountain environments on Earth (Barsch, 1996;Jones et al., 2018).
The long-term preservation of ice content in rock glaciers under permafrost conditions, situated under a rocky active layer, enable these landforms to creep over time periods of thousands of years (Haeberli et al., 1999;Krainer et al., 2015). The thermal conditions below the freezing point make the creep of ice-rich frozen ground inherently sensitive to climatic conditions and it is able to change over different time scales (Kääb et al., 2007;Delaloye et al., 2010;Sorg et al., 2015). Besides their 10 climatic and geomorphic significance, creeping frozen sediments can also be the source of rockfall and debris flows, and thus evolve into local-scale natural hazards (Kaufmann and Ladstädter, 2008;Roer et al., 2008;Stoffel and Huggel, 2012;Scapozza et al., 2014;Deline et al., 2015;Schoeneich et al., 2015;Bodin et al., 2016;Kummert et al., 2017;Scotti et al., 2017).
Creep speed of perennially frozen sediments is considered a valuable indicator of environmental, in particular of climatic and ground-thermal conditions. Creep speed of frozen debris in cold mountains is expected to increase with ground temperature 15 (Arenson and Springman, 2005a, b;Kääb et al., 2007;Arenson et al., 2015;Müller et al., 2016), although it is also influenced substantially by a number of other factors such as topography (e.g., slope), provision of ice and debris, temporal and spatial variations in ice content, rheology of the frozen debris, thickness of the creeping permafrost layer, as well as advection, infiltration, or internal production of water (Jansen and Hergarten, 2006;Cicoira et al., 2019a;Kenner et al., 2020). The role of water can enhance the response of creep speeds of frozen debris to rising ground temperatures when approaching melting 20 conditions (Kääb et al., 2007;Ikeda et al., 2008;Delaloye et al., 2010;Arenson et al., 2015;Hartl et al., 2016;Buchli et al., 2018;Cicoira et al., 2019b).
From the above, it becomes clear that inventorying contemporary motion and monitoring speed of bodies of ice-rich perennially frozen debris over time may contribute substantially to a better understanding of the climatic and geomorphic significance of such coherent creep of cold mountain slopes, its response to climatic changes, and the impacts from spatio-25 temporal variations of such creep, for instance changing slope and infrastructure stability, and alpine sediment transport.
Inventorying kinematics of rock glaciers and other bodies of frozen mountain sediments over regional scales seems currently best done using satellite radar interferometry, supported by optical and topographic data (IPA, 2020). This method enables detecting slow slope movement, not least rock glacier motion, over large regions and hundreds of individual landforms (Villarroel et al., 2018). The approach has been used to develop rock glacier motion inventories for the Swiss Alps (Strozzi et 30 al., 2004;Barboux et al., 2014), Sierra Nevada, California (Liu et al., 2013), Northern Iceland (Lilleoren et al., 2013), Brooks Range, Alaska (Rick et al., 2015), Carpathians, Romania (Necsoiu et al., 2016), north-eastern Tien Shan, Xinjiang/China (Wang et al., 2017), and the dry Andes of Argentina and Chile (Villarroel et al., 2018). These studies mainly relied on satellite L-band synthetic aperture radar (SAR) data such as from the JERS and PALSAR instruments, and C-band data such as from ERS and Sentinel-1. The categorisation of landforms into classified movement rates or orders of magnitude (cm/day, cm/month, dm/month, cm/year, etc.) has been proven useful (Barboux et al., 2014).
Rock glacier speed was shown to vary at a range of time scales, including seasonal, annual, decennial, or centennial (Delaloye et al., 2010;Scapozza et al., 2014). For a few rock glaciers worldwide, mostly concentrated in the European Alps, important 5 time series exist with seasonal to decennial time resolutions. Most of these measurements are based on terrestrial geodetic surveys. This implies that such time series cannot be established retrospectively for locations that have not been monitored in the past. Satellite radar interferometry is mainly available for 21 st -century creep velocities. The only means to measure 20 thcentury rock glacier surface velocities at yet unmeasured locations, and to reconstruct rock glacier speed variations at climatic time scales (i.e. at roughly decennial time scales) and back into the last century is by offset tracking in repeat historic and 10 modern air and satellite imagery (Kääb et al., 2007;Kaufmann, 2012;Sorg et al., 2015;Monnier and Kinnard, 2017).
The limited number of studies that have presented data on decennial variations in rock glacier speed in the past have drawn a contrasting picture of decadal-scale increases in speed, but also reported decreases or phases of decrease, as well as episodes of stable creep speeds (Schneider and Schneider, 2001;Janke, 2005;Krainer and He, 2006;Bodin et al., 2009;Delaloye et al., 2010;Kellerer-Pirklbauer and Kaufmann, 2012;Hartl et al., 2016;PERMOS, 2016;Eriksen et al., 2018). The long-term 15 observational studies underline qualitatively or quantitatively the influence of variations in mean annual air or ground temperature on the creep of perennially frozen debris, as suggested by theoretical and laboratory studies (see above). Decrease in speed was also shown to coincide with a significant surface lowering, thus hinting to a degradation of ice-rich permafrost (Bodin et al., 2009). Again, by far most of the above few studies on decennial variations in rock glacier speed stem from the European Alps. 20 The aim of the present study is to systematically inventory contemporary rock glacier creep speeds over two mountain ranges in the northern Tien Shan and to reconstruct speed changes since the mid-20 th century for selected landforms, with the specific goal to eventually provide spatio-temporal variations in the creep of frozen debris for a region outside the European Alps. The lack of availability of detailed long-term meteorological data distributed over the study region prevents a detailed analysis of the potential influence of temperatures, precipitation and snow cover on the observed rock glacier motions, though. 25

Study area and previous work
The mountain ranges of Ile Alatau and Kungöy Ala-Too (former names: Zailijskij and Kungej Alatau) are located at the border of Kyrgyzstan and Kazakhstan, Central Asia ( Fig. 1). They rise from the Kazakh steppe to the north and Lake Isyk-Köl (Issyk-Kul, 1608 m asl) to the south up to an elevation of almost 5000 m asl. Ongoing tectonic activity has caused many earthquakes with magnitudes of up to 8.0 and more (Yadav and Kulieshius, 1992;Lukk et al., 1995;Tibaldi and Graziotto, 1997), and the 30 stronger ones of these earthquakes have triggered a large number of mass movements (Delevaux et al., 2001;Korjenkov et al., 2004;Strom and Korup, 2006). The dominant lithology of the mountain ridges is granite, but weaker schists build the lowerlying saddles as well as larger parts of the valleys.
The climate of the study region is continental with mean annual precipitation exceeding 1000 mm on the northern slopes, but less than 800 mm on the southern ridges and in general for elevations above 3000 m asl (Bolch, 2007;Narama et al., 2010).
Minimum precipitation is typically recorded during winter under the Siberian High. Maximum precipitation is found for spring 5 and early summer. The zero-degree isotherm is estimated at slightly above 2700 m asl. (Bolch, 2007). Ground temperature measurements showed that permafrost occurrence is likely above 3200 m asl. (Gorbunov, 1996). Air temperatures have been increasing since the 1970s with annual rates of about 0.01 °C, and higher rates during late summer and winter (Bolch, 2007;Marchenko et al., 2007;Narama et al., 2010;Sorg et al., 2012). Ground temperature measurements at depths of 5-25 m at about 3300 m asl. show an increase between 0.3 °C and 0.6 °C from the 1970s to 2004 (Marchenko et al., 2007). All rock 10 glaciers investigated here in more detail using photogrammetry (Gorodetzky, Morenny, Archaly, Ordzhonikidze, Karakoram and Kugalan Tash) originate in areas where permafrost is very likely whereas the fronts of the large rock glaciers are often located at elevations slightly below the lower boundary of regional permafrost occurrence (Bolch and Gorbunov, 2014). In the hemispheric-scale 1km-resolution permafrost model by Obu et al. (2019), all the rock glaciers studied in detail have permafrost probabilities > 0.5, and modelled mean annual ground temperatures (MAGT) between -0 and -2 °C. In line with the results 15 from Bolch and Gorbunov (2014), the model indicates that Karakoram rock glacier reaches down into non-permafrost terrain (MAGT ~+1.4 °C).
Rock glaciers are widespread in the study area, and some are remarkably large with areas exceeding 2 km². A first rock glacier inventory was compiled based on aerial photography from the 1970s. More than 1000 rock glaciers (amongst them about 180 inactive ones according to visual interpretation by the below authors) were identified in Ile Alatau and Kungöy Ala-Too, 20 covering a total area of about 90 km² (Gorbunov and Titkov, 1989;Gorbunov et al., 1998). The average altitude of the active rock glaciers is 3400 m asl., the lowest elevation of an active rock glacier front is less than 2700 m asl. Most rock glaciers in Ile Alatau are on north-facing slopes, whereas in Kungöy Ala Too the southern aspect prevails. The highest rock glacier density is found in the central and western part of Ile Alatau as well as in some north-facing side branches of the Chon-Kemin and Chilik valleys where the percentage of glaciated terrain is 10-20% (Titkov, 1988;Gorbunov and Titkov, 1989). Rock glacier 25 density decreases towards the outer parts of the mountain ranges due to decreasing elevations of the ridges. In a more recent inventory, Bolch and Gorbunov (2014) identified 72 rock glaciers in the central part of Ile and Kungöy Ala Too based on Landsat ETM+ imagery and aided by aerial images as well as high-resolution satellite images from GoogleEarth. Blöthe et al.
(2019) mapped rock glaciers over High Mountain Asia, including our study region, with the purpose to statistically model and detect those landforms that divert, block or else impact on river channels. One of our study rock glaciers, Ordzhonikidze, is 30 explicitly mentioned in their study, but also Karakoram and Kugalan Tash rock glaciers studied here clearly divert the river channels at their snout.
Displacement measurements of large boulders on a few large rock glaciers using ground surveys and aerial imagery from 1969, 1979 and 1994 revealed annual speeds between less than 1m/a to extraordinarily high rates of more than 11 m/a (Gorbunov et al., 1992). On average, there was an increase of the surface speed found between the periods 1969 -1979 and 1979 -1994, except for the Karakorum rock glacier, which is at the lowest elevation. Displacement rates of individual boulders vary greatly within rock glaciers. There are some areas where the latest measurements show inactivity (Gorbunov et al., 1992;Kokarev et al., 1997). Most of the measured boulders showed displacements in slope direction, but those rock glaciers reaching the valley floor also exhibit lateral displacement, pointing to a spread onto the valley floor. It was suggested that some rock glaciers in 5 Tien Shan are subject to rapid advances, maybe caused by earthquakes. This might be the reason why some rock glaciers have been shown to advance into forested areas (Gorbunov, 1983). Four of the rock glaciers studied in Gorbunov et al. (1992) are also investigated in detail in this study, and results are compared in the discussion section of this contribution. Sorg et al. (2015) used tree-ring growth anomalies from 280 trees growing on rock glaciers in the study area, or close to it, as an indicator for rock glacier activity and compared their results to climatic data and detailed glacier mass balances for the 10 period 1895 -2011. They find that the activity of the rock glaciers studied correlates well with decadal variations of summer air temperatures and glacier mass balancerock glacier activity/tree growth anomalies are positively linked to high summer temperatures and inversely to glacier mass balances. As a secondary focus, these authors also use preliminary offset tracking between repeat aerial and satellite imagery to show that on a multi-decadal scale the variations in rock glacier speeds show more contrasting and site-specific changes, with accelerations and decelerations, both varying over time and even within the 15 same rock glacier system. They attribute this behaviour partially to the different stages of the rock glaciers on their path towards inactivation under continued atmospheric warming in the region. Sorg et al. (2015) did not use radar interferometry in their study nor compile a kinematic rock glacier inventory. It is the offset tracking work started in Sorg et al. (2015) that we extend in the second part of the present study in terms of data and time steps, number of rock glaciers, and methods used. The specific rock glaciers studied to this end are Morenny, Archaly, Gorodetzky, Ordzhonikidze, Kugalan Tash and Karakorum (Fig. 1). 20 These rock glaciers were chosen because they were the focus of previous studies including rock glacier monitoring work, and they were the only ones with suitable airphotos for time-series measurement available to us. More details about these rock glaciers are introduced under the results section.

Radar interferometry 25
An inventory of rock glaciers and other periglacial processes in the Ile and Kungöy Ranges of Northern Tien Shan was in this study compiled by visual analysis of differential SAR interferograms and satellite optical images. We computed differential SAR interferograms (DInSAR) with short baselines and time intervals between 1 day and one year from the ERS-1/2 tandem mission (1998)(1999) Table S1). Images acquired along both ascending and descending geometries and during summer (snow-free) 30 and winter (dry snow) conditions were employed. For topographic reference and orthorectification, we computed a DEM from TanDEM-X data at a spatial resolution of 10 m. Bistatic TanDEM-X acquisitions along ascending (13 Jan 2012 and 24 Jan 2012) and descending (15 Dec 2013 and6 Jan 2014) orbits were combined to reduce problems in areas masked by layover and shadow, which differently affect opposite slopes along north-south oriented valleys in mountainous regions (Ambrosi et al., 2018). In general, we estimate the errors due to atmospheric distortions and noise in the DInSAR images to be in the order of 1/4 of a phase cycle (i.e. 1/8 of the radar wavelength used) (Strozzi et al., 2001). All processing was done using the GAMMA radar software. 5 By combining interpretation of all interferograms, signals related to the motion of periglacial phenomena over various time intervals were first identified from the wrapped interferometric signals following the recommendations given in Barboux et al. (2014) and marked as point features. In a second step, the extent of the InSAR movement signals was mapped based on the interferograms, and in addition supported by optical high-resolution imagery from GoogleEarth/Map and Bing Maps. Next, the deformation rates of the mapped landforms were classified using five classes (i.e. 0-2, 2-10, 10-50, 50-100 and >100 cm/a). 10 In short, interferogram interpretation and estimate of movement class is based on a large stack of interferograms, not only 12day Sentinel-1 interferograms as shown in the figures. In suitable interferograms, movement fringes are counted for moving areas with phase coherence and this fringe count is converted to line-of-sight movement using the respective radar wavelength.
Confined areas of coherence loss (fringe noise) can be due to high surface speeds that decorrelate the interferometric signal.
The choice of classes is routed in guidelines from the Swiss government for landslide, rockfall and debris flow hazard 15 management, which was the only such recommendation with respect to radar interferometry available at the time of inventory compilation, and is certainly open to discussion. An updated version of the openly available inventory will rather follow classes meanwhile suggested by an international working group (IPA, 2020). In cases where the estimated velocity was close to the upper boundary of a velocity class, the polygon was assigned to this faster class as the 1-dimensional line-of-sight measurement provided by radar interferometry represents only one movement component and thus typically underestimates the actual three-20 dimensional surface movement. The same was done with natural temporal variations of surface displacement rates, i.e. if two or more classes were present during the observation time-span, the highest displacement rate was used to determine the velocity class. Because the interferometry-derived displacement rates are measured along the satellite line-of-sight, we empirically considered the line of sight in the assessment of the state of activity of the mapped landforms (see above). Approaches to rock glacier kinematics inventorying from radar interferometry are currently developing much (IPA, 2020). At the time of 25 compilation of our inventory, it was practice to rely on expert decisions based on wrapped interferograms. For current compilation of such inventories, recent computational and data availability advances are exploited to further process interferograms, for instance accounting for variable line-of-sights or exact time intervals.
In the final step of our analysis, the specific interferometric fringe pattern, optical images and the shaded relief of the DEM were visually interpreted to classify the process types associated to the observed slope movements. In addition to the rock 30 glaciers as main objective of this study, movements related to solifluction / debris movements, dead ice / subsidence, debriscovered glaciers and ice-cored moraines were classified. Related recommendations from an international working group (IPA, 2020) were not yet available at the time of inventory compilation. High-resolution satellite images were the main source for the classification into the above classes. In cases of missing information due to low resolution or cloud cover in the images, the interferograms were consulted as secondary source for the interpretation. The following exemplifies some considerations taken during the expert decision. For instance, rock glaciers usually show an overall increase of speed towards the front and inner region of the polygon mapped, which helps to identify them against other movements. Rock glaciers can be discriminated from debris-covered glaciers by thermokarstic features such as supraglacial ponds, which lead to loss of radar phase coherence and are visible in optical or even radar data. In addition, ice is sometimes visible in high-resolution optical images. Melting 5 dead ice bodies have characteristic spot-like fringe patterns from local subsidence that are visible only in late summer interferograms. Similar appearance in the interferograms can sometimes be observed for debris movements (soli-/gelifluction).
This process type is, though, typically found on inclined slopes but lacks large-scale compressional features (lobes). Ice-cored moraines are connected to the former lateral or terminal moraines of a retreated glacier, a situation that can often be recognized in all kinds of image and topographic data. For cases where the InSAR signal could not be unambiguously assigned to a process 10 type and where also optical images and DEMs were not of help or not available with sufficient quality, an "undefined" class was introduced.
Above the maximum inventory movement rate of 100 cm/a, motion could not be quantified anymore, but detection was still well possible. This maximum measurable displacement rate does not only depend on the wavelength used, but also on the movement type (e.g., spatially coherent vs. locally variable) and terrain (e.g., terrain slope relative to radar incidence angle, or 15 movement azimuth relative to satellite azimuth). In the case of Sentinel-1 interferograms, for instance, assuming a maximum measurable displacement rate of half a wavelength (2.8 cm) (Villarroel et al., 2018) gives for the shortest temporal baselines available over the study area (12 days) a maximum value of 85 cm/a. One reason for classifying decorrelated areas in the > 100 cm/a class is that motion decorrelation usually happens actually under higher displacements as this minimum value of 85 cm/a. A second reason is that the topographic correction factor is in realistic cases not lower than 1.1-1.2, making a line-of-20 sight displacement rate of 85 cm/a rather belonging to a landform of > 100 cm/a class. Data from the ERS 1-day tandem phase were not good enough, likely due to snow melt and snow changes, to quantify motion in areas that were decorrelated in interferograms from other sensors.
Currently, methods for computing seasonal time series of rock glacier speed from repeat radar interferometry are being developed (Strozzi et al., 2020), but these methods and time series of suitable data were not yet fully available at compilation 25 of the kinematic inventory of this study. Further, producing such time series for a number of rock glaciers would have by far exceeded the resources available for our study.

Matching of repeat optical images
For four sites with rock glaciers in the study area (Gorodetsky/Morenny/Archaly; Ordzhonokidze; Kugalan Tash; Karakorum; radiometric quality, distortions and sensor types used in this study are very different. Roughly, we used three categories of data: historic airphotos, historic spy satellite images, and contemporary high-resolution satellite data (Supplementary Table   S2). The contemporary satellite images from the DigitalGlobe (now Maxar) constellation or Pléiades were oriented using standard procedures based on the provided rational polynomial function coefficients, and checked with a few control points extracted from GoogleEarth. Over the Gorodetsky/Morenny/Archaly rock glacier group, Pléiades stereo data were available.
A corresponding digital elevation model (DEM) was therefore constructed, and the satellite images were orthorectified based on this DEM. For the other three sites without own DEMs produced, orthoimages were computed using the C-band SRTM DEM (30-m resolution), the ALOS PRISM AW3D DEM (commercial version, 5-m resolution version), and the TanDEM-X  5 WorldDEM (commercial version, 12-m resolution version). The High-Mountain-Asia DEM (Shean, 2017)produced with stereo data from the DigitalGlobe/Maxar archivewas also investigated; it turned out to be very good for the Karakoram and Ordzhonikidze rock glaciers, but showed large voids over the Kugalan Tash rock glacier. Based on the analysis of these different DEMs such as from DEM differences, level of detail, and void cover, and based on the distortion patterns over stable ground between the corresponding orthoimages produced, it was decided to use the TanDEM-X WorldDEM for the generation 10 of the orthoimages, except for Gorodetsky/Morenny/Archaly where the Pléiades DEM was used. At the time of processing, the TanDEM-X WorldDEM appeared to be the most complete, detailed and consistent DEM over all sites. The orthoimages generated from the contemporary optical satellite data served as source of georeference for the older images.
The airphotos used in this study stem from Soviet-era archives, and the digital copies available to us were of very variable quality, typically scans of copies of the original images using consumer-grade scanners, or similar. The images had thus 15 significant distortions and very variable radiometric contrast. No camera calibrations were available. In addition, we used highresolution declassified spy satellite images from Corona missions (KH-4A and KH-4B) (Altmaier and Kany, 2002;Bolch et al., 2008;Narama et al., 2010;Goerlich et al., 2017). These images have a complex panoramic geometry, in addition to film distortions. Thus, to both the airphotos and Corona data, no rigorous sensor models were applied but rather empirical models were used consisting of 4 th -order rational polynomials computed from around 40-50 ground control points at each site. Ground 20 control was extracted from the contemporary satellite orthoimages and DEMs. All processing of the optical data was done using the PCI Geomatica software (version 2018sp1).
Eventually, horizontal displacements were measured between pairs of the orthoimages using standard subpixel-level normalized cross-correlation (Kääb and Vollmer, 2000;Debella-Gilo and Kääb, 2011;Kääb, 2014). For the contemporary satellite images and some particularly good 20 th -century images, entire flow fields could be measured successfully. For the 25 majority of the older image data, the resulting velocity fields contained substantial areas were image matching failed, as indicated by low correlation values, or vectors pointing randomly at different directions. Time series were only constructed for points that had reliable measurements for all periods, i.e. measurements that passed thresholds on correlation coefficients (depending on individual image quality, but usually between 0.4 and 0.6) and manual removal of obvious outliers. The time series shown in the following results and figures are average time series for all individual measurements within manually 30 delineated zones. These zones are also shown in the figures. The purpose of creating such time-series zones was to focus on central rock glacier parts with substantial movement, and to exclude rock glacier parts that were considered less representative, for instance due to location at the margins of the rock glacier (system) or due to uncertain measurements. Even for these averaging zones, the low number of suitable measurements in the historical images typically limited the choice of clusters of points with available time series. The accuracy of temporal speed variations was then evaluated using three measures: (i) speeds on stable ground (for Gorodetsky, Morenny, Archaly, and recent periods for the other rock glaciers only), (ii) variation of speeds within the point clusters (spatial variance), and (iii) variation of speed variations between periods and clusters (spatiotemporal variance). Potential uncertainties in the time series are thus given as the normalized median absolute deviation (NMAD; Höhle and Höhle, 2009) of stable ground speeds (red error bars, (i)), as the standard deviation of speeds over one 5 period (blue error bars, (ii)), and as standard deviation of speed differences between two periods (grey error bars, (iii)). The purpose of error bars (ii) is to visualize the homogeneity of speeds within a cluster. Variations of this measure over a time series give then hints to measurement uncertainties, under the assumption that the velocity pattern within the cluster is not changing strongly over time. The error budget resulting from the highly variable quality of the base images, their georeferencing, from the DEMs used for their orthorectification, and the offset tracking procedure is examined in more detail 10 in the discussion section of this study. The reason why stable ground offsets were included as error bars for Gorodetsky, Morenny and Archaly rock glaciers only is that these three landforms are surrounded by large areas of flat terrain, a surface type and topography that we considered representative for the rock glaciers themselves. The other rock glaciers are in large parts surrounded by steep terrain with a large number of image matching problems and larger impacts from DEM errors on the image orientation and orthorectification, in particular for the airphoto and Corona data. Automatic processing of robust 15 and meaningful stable ground error bars was thus not feasible for the time series of the other rock glaciers, and respective errors are given in the result text where available.
In a supplementary analysis, we also relate the advance rates of rock glaciers to the surface speeds above the fronts. The ratio between both quantities is a function of the ice content that melts out at the front, and the vertical profile of horizontal speeds, i.e. the total material flux into the front (Kääb and Reichmuth, 2005). We roughly estimate the advance rates by measuring 20 some ten manually selected distances along the main flow direction between repeat positions of the front base using flickering orthoimages of modern-era high-resolution satellite images. It is clear that such estimate of advance is coarse, and difficult anyway due to local variability from e.g. falling boulders and erosions. Our aim is thus only to investigate large differences between rock glaciers in a qualitative way. Detailed quantitative analysis of such advance rates would require disentangling the compound effects of ice content and vertical velocity profile, which is to our best knowledge not possible using air and 25 satellite images. Respective conclusions are only possible when the ratio between advance and surface speed is very small (approaching zero; pointing to a very high ice content melting out and/or very shallow moving layer, both preventing strong advance), or approaching 1 (very low ice content and little decrease of horizontal speeds with depth, i.e. so-called block movement). Also, the resolution and quality of the image data available to us is not sufficient to measure time series of advance rates. 30 Extending photogrammetric time series to other rock glaciers in the region would certainly be of value, but the large difficulties of acquiring Soviet era airphotos, the high costs of the commercial modern high resolution satellite images, and the large effort necessary to process the data consistently from very different sensors and with very different qualities, prevented to extend the time series to even more rock glaciers as already done in this study.

Rock glacier motion inventory
The inventory of rock glacier creep speed and movement of periglacial and glacial processes over the Ile and Kungöy Ranges of the northern Tien Shan covers a surface of more than 3000 km 2 and includes more than 900 landforms, namely 551 rock glaciers, 184 solifluction features / debris movements, 46 dead-ice / subsidence forms, 18 debris-covered glaciers, 14 ice-cored 5 moraines (for inventory access see data availability section). Finally, 98 features remain undefined due to the lack of good quality high-resolution images available at the time the inventory was created. Of the active rock glaciers, which are homogeneously distributed over the study region, 45 are characterised by a rate of motion exceeding 100 cm/a. We identified 21 large rock glaciers exceeding 1 km 2 in size with the largest being about 2 km 2 . The mean rock glacier size is 0.27 km 2 .
Some of the largest rock glaciers have also the largest elevation difference of up to 1,500 m. Several rock glaciers extend 10 below 3000 m asl., with the lowermost elevation of active rock glacier fronts being at about 2700 m asl., amongst them Karakorum and Ordzhonikidze rock glaciers, which are analysed in more detail in this study. The highest rooting zone of a rock glacier is situated above 4000 m and the mean elevation is 3480 m. Six rock glaciers within the class > 100 cm/a are further analysed in this study by matching of repeat optical images: Gorodetsky/Morenny/Archaly, Ordzhonikidze, Karakoram, and Kugalan Tash (Fig. 1). 15 A selected Sentinel-1 interferogram for all these rock glaciers demonstrates how the motion cannot be quantified over large parts of the rock glaciers as no consistent movement fringes are present (Fig. 2). However, the moving areas can still be detected and mapped even if movement fringes cannot be reliably unwrapped (i.e. integrated) into displacement rates or coherence is completely lost.

Gorodetsky, Morenny and Archaly
Surface velocities from matching repeat images, in the following termed in short "photogrammetric velocities", between 2012 and 2016 (for exact dates and image sources see Table S2) over Gorodetsky rock glacier (Fig. 3a) show a pattern coherent 10 viscous flow, with several zones of comparably high speed (Fig. 3c). Highest speeds of up to 1.4 m/a are found along the lower southern margin of the rock glacier, which agrees well with the interferograms (Fig. 2) where strongest signs of coherence loss (noisy fringes) are found over the same zone. The NMAD of stable ground offsets for 2012-2016 is 0.25 m (0.06 m/a).
Gorodetsky rock glacier seems to be connected to the forefields of five glaciers, three of which are visible in Figure 3a, whereas two are located more to the east. The zones between the present-day glaciers and the main rock glacier do not show any 15 significant movement. Close inspection of the repeat high-resolution optical images and the interferogram Fig. 2b rather reveals a mixture of areas characterized by rapidly changing thermokarst features, such as melt ponds, melt depressions, or small slope instabilities, and areas that appear stable and without any visible, ice-related activity. Together, these lines of evidence suggest that the rock glacier has today no direct contact with the glaciers.
A time series of speeds for a cluster of points (ca. 45) in the central portion of the rock glacier (white dashed outline in Fig.  5 3c) suggests a decrease of speed over 1953-1964-1971-2012. By contrast, over the period 2012-2016-2018, we measured a slight increase in speed, though hardly exceeding the error bars. The large gap of data between 1971 and 2012 prevents a more detailed assessment and only enables measurement of average speeds over a large part of the observation period.
The results for Archaly and Morenny rock glaciers are comparable to those obtained for Gorodetsky, except that the 1953 airphoto covers only Gorodetsky. The forefields connecting the rock glaciers and the glaciers above them (Fig. 4c) show clear 10 signs of thermokarstic ice degradation and slow/no movement as well, but the evidence is spatially less extensive than at Gorodetsky. A distinct fast-flowing branch (up to 2.5 m/a) diverts from the mainly south-north flowing Morenny rock glacier towards the northwest (Fig. 4e). The InSAR results (Fig. 2) are consistent with the velocity field ( Fig. 4e) in that they show coherent fringes over much of the slower-moving part of the rock glacier, but strong noise from phase coherence loss over the fast-moving lateral branch. Overall, Archaly rock glacier moves faster than Morenny (Fig. 4d), which is also reflected in the 15 more decorrelated interferometric phases (Fig. 2). Though, also on the central part of Morenny rock glacier speeds are in the order of >100 cm/a, in a zone just above the time series cluster, i.e. consistent with the interferogram interpretations (Fig. 1).
The time series for Archaly (~40 points) and Morenny (~75 points) rock glaciers ( Fig. 4a and b) look similar to Gorodetsky for 1964 onwards, though with larger speed increases relative to their error bars.

Karakoram
The velocity field on Karakoram rock glacier (Fig. 5a,c) shows one coherent stream of material from a partially ice-covered rock wall (Fig. 5a). A few thermokarst depressions and lakes are found in the root zone of the rock glacier. The lowermost part of the rock glacier hosts significant vegetation including dozens of trees (Sorg et al., 2015), visible in red colour in the infrared 5 false colour composite Fig. 5a. Maximum speeds recorded over 2009-2017 reach up to 6.5 m/a. Speeds in the upper part of the rock glacier, in the zone of thermokarst features, are in the order of a few dm/a. The NMAD of stable ground offsets for 2009-2017 is 0.65 m (0.08 m/a), but matching quality around the rock glacier is clearly lower than on the rock glacier. The high speeds of Karakoram rock glacier are consistent with the interferogram given in Fig. 2, where strong phase noise is visible.
The high speeds and associated surface changes result in a situation in which a larger number of points with valid velocity 10 measurements over all measurement periods can only be found in the lowermost portion of the rock glacier. There, the average of the point cluster (ca. 20 points) shows a steadily increasing speed over the period 1960-1971-1980-1988-2001-2009-2017 ( Fig. 5b). Available airphotos from 1977, 1985 and 2008 have not been used in this time series because they were of too low quality to enable reliable measurements, or their dates were too close to other image dates so that measurements did not exceed the measurement noise level. The advance rate of Karakoram rock glacier can only be measured in a meaningful way in the 15 north-western portion of the rock glacier front as the other parts of the front terminate in a river. In this north-western part of the front, the advance rate (~1 m/a) is very similar to the surface speeds above the front, leading to a ratio close to 1. This points to block movement of the entire rock glacier column (i.e. little decrease in horizontal speeds with depth until a basal shear horizon) and very little ice content left at the front.  Table S2). Background image: GeoEye, 29 Jul 2017, © Maxar.

Ordzhonikidze 10
Ordzhonikidze rock glacier is a system of several streams of viscous flow, all connected to small cirque glaciers or their forefields (Fig. 6a). The lowermost part of the longest stream reaches into forested terrain and carries vegetation, including dozens of trees (showing in red in Fig. 6a) (Sorg et al., 2015). For all individual rock glacier streams, a coherent flux of material seems to exist from the cirque glaciers into the rock glaciers, with surface speeds of around 1 m/a. Highest speeds of up to 4 m/a over 2016-2018 are found on a steep narrow part of the longest stream. The NMAD of stable ground offsets 2016-2018 is 0.65 m (0.32 m/a), but matching quality around the rock glacier is clearly lower than on the rock glacier. The interferometric fringe pattern (Fig. 2) and the interferometric interpretations ( Fig. 1) are consistent with photogrammetric measurements (Fig.   6e), with phase noise over the fastest rock glacier parts, and more coherent fringes over the slower parts. Velocity time series 5 were computed for four zones (A: ~ 60 points; B: ~100; C: ~70; D: ~50; Fig. 6e) and over the periods 1953-1966-1971-1985-2009-2013-2016-2018 (for exact image dates and types see Tab. S2). A period of higher speeds over 1966-1971 is found for all four zones (Fig. 6b-d). Although with variable timing and to variable degrees, we detect an overall increase in speed for all rock glacier zones investigated. Rock glacier advance rate ratios are difficult to characterize because of the many separated fronts. Roughly, at the front below zone A, advance rates are in the order of 0.7 m/a and speeds above this front are between 10 1 and 2 m/a (corresponding to advance ratios of around 0.35 -0.7). At the front below zone D, the advance is in the order of 0.7 m/a, and surface speeds are around 1.2 m/a (thus a ratio of around 0.6), based on advance rates and speeds over the period 2009-2018.

Kugalan Tash
Kugalan Tash rock glacier is connected to a cirque with no significant ice bodies remaining at the present day (Fig. 7a). Our earliest airphotos of 1956 show clearly more ice cover than today, both in the rock walls and at their base. Movement rates between this cirque and the rock glacier are smaller than the error bars. Only very little vegetation is growing on the rock glacier ( Fig. 7a) (Sorg et al., 2015). Maximum speeds found on the stream of viscous flow over 2012-2014 are around 2.5 m/a. 15 The NMAD of stable ground offsets 2012-2014 is 0.65 m, like for Karakoram and Ordzonikidze (0.32 m/a). The photogrammetrically measured rock glacier stream is well reflected in the interferogram depicted in Fig. 2, with coherent fringes in the slower parts, and phase noise due to coherence loss over the faster parts. Speed time series for 1956-1973-1988-2004-2012-2014 were measured for two zones (A: ~20 points; B: ~30). Zone A shows a steady increase in speeds (Fig. 7c).

Radar interferometry versus photogrammetry
The photogrammetric velocity fields derived from repeat airphotos and high-resolution optical satellite data in this study agree reasonably well with the radar interferograms and their interpretation. Noteworthy, both datasets have been obtained 15 independently as the photogrammetric velocity fields were not accessible to the team processing and interpreting radar interferometry data, and vice versa. Nonetheless, all rock glaciers investigated using photogrammetry were classified correctly from the interferograms (Fig. 1). When having a closer look at the full resolution of the radar interferograms, however, small differences can be seen. As radar interferometry provides in principle only one motion component (line-of-sight)by contrast to photogrammetry where two horizontal components are obtainedany discrimination of lateral from vertical movements is 5 obviously more difficult in radar interferograms. By way of example, a small elongated zone of phase coherence loss (phase noise) can be detected in the eastern sector of the main part of Gorodetsky rock glacier (Fig. 2), and has been interpreted as a fast-moving zone within the rock glacier (Fig. 1). Photogrammetric investigation, however, shows that this zone should be interpreted rather as a zone of thermokarstic erosion, a process that can cause loss of interferometric phase coherence. One needs to underline, however, that radar interferometry and multitemporal photogrammetry require very different datasets and 10 that the two approaches usually serve different purposes. Radar interferometry is typically used to reliably map the activity of a large number of rock glaciers over large areas, a target that can only be achieved by using free, weather-independent, and easily accessible baseline data (e.g. Sentinel-1). By contrast, photogrammetric approaches as the ones used here, require very specific, highly-resolved data that cover only small areas, and can be costly or else difficult to obtain.
In addition, we would also like to stress that the motion patterns seen in the wrapped radar interferograms (Fig. 2) agree very 15 well with the photogrammetric velocity fields. Fast moving rock glacier parts are characterized by interferometric phase noise, whereas slower-moving parts show more coherent interferometric fringes (e.g., Fig. 8, and Supplement for the other rock glaciers). In the 12-day interferograms displayed in Figs. 2 and 8, the transition between both cases seems to be indeed at speeds of around 85-100 cm/a, which is in good agreement with estimates from Villarroel et al. (2018).
Whereas our motion classes are based on quantitative measurements, it is clear that our landform classification involves some 20 uncertainty especially for process types that naturally have no sharp boundary between them (in spatial terms or in terms of process definition) but rather transitions, such as the downstream transition from a debris-covered glacier into a rock glacier.
Further it is clear, that the radar and optical data available at the time of classification, for instance their resolution, will influence the class decision.

Uncertainties in the photogrammetric measurements
In this section, the error budget of our photogrammetric displacement measurements, and thus of our speed time series, is analysed. According error budgets are typically composed of (a) overall shifts between the orthorectified data, (b) lateral shifts 10 in the orthoimages due to errors in the DEM used for orthorectification, (c) distortions in the raw data or in the sensor model that propagate into the orthoimages, and (d) image matching uncertainties and errors (Kääb et al., 2016;Kääb et al., 2017Kääb et al., , 2019. Biases of type (a) are much minimized by our processing scheme and should not influence our results in a significant way.
One orthorectified high resolution satellite scene is used as master for georeferencing of all other images over a site, and image 15 sections are additionally co-registered over stable ground during the offset tracking procedure.
Orthorectification offsets due to DEM errors, i.e. error budget component (b), are difficult to quantify precisely. They will in particular affect the airphotos we used as those have typically much larger opening angles (field of view versus flying height) compared to the satellite images. In the very worst case, a DEM error could translate 1:1 to a lateral offset for off-nadir projection angles of 45°, which are typical at the margins of airphotos. We can estimate such orthorectification offsets from stable ground offsets (e.g. Figs. 3b and 4a-b), but we cannot be sure that these estimates are representative for the rock glacier 5 surfaces as the elevations there could have changed over time with respect to the DEM used for orthorectification (here the TanDEM-X WorldDEM), or DEM errors on rock glaciers could be systematically different from DEM errors outside the rock glaciers. In this study, we computed a number of orthoimages using different DEMs (see Section 3.2) and analysed offsets in them in order to select the best DEM. From the stable ground offset tests we performed, we find NMADs of offsets around the rock glaciers ranging from around 1-3 m for the old airphotos and Corona data, and 0.25-0.65 m for the modern satellite 10 images. Another, more stringent way to estimate orthorectification errors (b) would be to analyse offsets between stereo orthoimages, i.e. different images of the same scene (typically stereo images) rectified by the same DEM (Kääb and Vollmer, 2000;Kääb et al., 2016;Altena and Kääb, 2017). We had no suitable image data available to perform such tests, though.
Distortions in the raw data or sensor model used, i.e. errors of type (c), are especially expected for the Corona spy satellite data and some of the low-quality airphotos obtained. We have attempted to minimize such effects by using 4 th -order rational 15 polynomial sensor models computed from around 40-50 ground control points per image. These errors can be estimated from stable ground offsets as well, and are also contained in the respective numbers given above. By contrast to the orthorectification offsets, it can be assumed that image or model distortions are independent of the image content, and affect thus stable areas in a similar way than the rock glaciers.
The final error category (d) concerns matching errors. These can range from gross outliers to fine matching inaccuracies. Gross 20 outliers in our displacement measurements are eliminated by thresholds on the correlation coefficient (typically between 0.4 and 0.6) and a manual selection of points that have consistent measurements in all measurements periods for a site. Thereby, those vectors appearing inconsistent with respect to surrounding vectors (i.e. mainly sections of the vector field with random variation of directions) have been disregarded. In addition, we made sure that the time series obtained from mean speeds for point clusters were consistent with the time series obtained from median values for the point clusters, as the latter are less 25 sensitive to outliers. Mean values give better representation of the actual measurements, but the median values are a suitable check if outliers were properly removed by the thresholds and visual removals described in section 3.2. The remaining matching accuracy varies certainly much with the corresponding image quality, including resolution, contrast, and/or noise.
Effects from matching inaccuracies are minimized by averaging speeds of point clusters over several tens of points. Again, this matching accuracy can be estimated from stable ground offsets, in particular for surfaces that have similar visual contrast 30 as the rock glaciers, and is included in the numbers given above under error component (b). A further possibility for describing the uncertainty of individual average speeds for the time series intervals could be average correlation values for each step.
Correlation values alone are, however, not a very unique measure to characterize (geometric) matching uncertainty. Rather they indicate radiometric similarity and signal outliers, which we, however, detected visually in this work.
In stable ground offsets, the combination of error types (b)-(d) is present. For all years for Gorodetsky, Morenny and Archaly rock glaciers, and one recent interval for the other rock glaciers, we compared the normalized median absolute deviations (NMAD) of stable ground offsets to the standard deviations of speeds within point clusters, and to the standard deviation of speed differences to find that the latter two seem to be on overall a more conservative uncertainty measure, which we thus show for the other intervals of the time series. The latter two standard deviations contain also natural and, thus, non-erroneous 5 variability in speeds and speed differences. On the other hand, we know that orthorectification errors (b) might be larger on the rock glaciers than for flat areas around them, so that we prefer to use an error measure that is more conservative than the stable ground offsets.

Speed time series
The photogrammetric velocity fields from this study show a wide spatial variation of rock glacier speeds, with the maximum 10 speeds of up to 6.5 m/a found on Karakoram rock glacier. Whereas the Karakoram, Ordzhonikidze, and Kugalan Tash rock glaciers, or rock glacier systems, consist of well identifiable individual streams, Gorodetsky, Morenny and Archaly rock glaciers display more complex flow patterns, viscously spreading out over a flat valley bottom.
This study also shows that the temporal variations of rock glacier speeds from repeat airphotos and satellite data, available since the 1950s, displays a contrasting picture. The only feature in the time series common to all rock glaciers studied in more 15 detail is indeed an increase in speeds in recent years, but results clearly differ with respect to the ratio between the amount of speed increases and their error bars.
Karakoram rock glacier shows a consistent increase in speed over the entire study period, but due to its fast movement, which much limits reliable measurements over all measurement periods and all sections, this can only be confirmed for the lowermost part of the rock glacier. Through the manual tracking of 15 individual boulders in airphotos over the period 1969period -1977period -1984period , 20 Gorbunov et al. (1992 find a similar behaviour of speed increase in the lowermost part of the rock glacier (their measurement points 1 and 3) (Gorbunov and Titkov, 1989). For the parts further up, these authors documented a decrease in speed over their two study periods, interpreting these changes as a wave-like acceleration/deceleration process as the rock glacier material is travelling downstream. Measurements on our images over these rock glacier parts are too scarce and have too large error bars to allow for a meaningful comparison. The lowermost rock glacier part, where our measurement points are situated, deforms 25 under a regime of strong longitudinal compression, with longitudinal strain rates in the order of 0.012/a derived from the photogrammetric velocity fields, and consistent with the many transverse ridges visible (Kääb and Weber, 2004). We suggest this rock glacier part responds to variations in mass flux from above through a combination of variations in horizontal speeds and surface elevation changes. Under the assumption of incompressible material these two components of rock glacier kinematics would be correlated, but we have no measurements of elevation changes available to disentangle them. As indicated 30 by the advance ratio of close to 1 (advance rate and surface speed both in the order of ~1 m/a), the ice content of the front could be quite low. Both for the lowermost rock glacier part (Fig. 5a) and the few points measurable in the middle part, speeds have more than doubled over 1964-2017. The results from the present study are also consistent with those from Sorg et al. (2015), with both studies showing a steady increase in speed over time. The measurement points used in Sorg et al. (2015) extend not as far down to lower rock glacier elevations as the increased number of points in the present study, so that the speeds found in their work average to higher values compared to the present study.
Kugalan Tash rock glacier shows an overall increase of speeds over 1956-2014 as well, but at a lower significance level than for Karakoram rock glacier (Fig. 7c). The increase in speeds on the lowermost rock glacier part (zone B) is lower than for zone 5 A. Zone B, under a strong compressive creep regime (up to about 0.03/a), reaches and diverts the river, similar to the situation observed for the lowermost part of Karakoram rock glacier. The photogrammetric measurements in Sorg et al. (2015) show, by contrast to the reprocessing and extension performed here, higher speeds during 1956-1971, a feature that we cannot confirm in the new assessment anymore. Through use of an improved DEM in our new processing and other improvements related to co-registration of images we consider the results presented here to be more reliable. 10 Overall, all investigated zones on the Ordzhonikidze rock glacier system show an acceleration over the period 1953-2018 as well, a phenomenon that is most clearly visible in the higher-lying northern rock glacier (zone A in Fig. 6e). A striking peak in speeds is visible over the period 1966-1971, for which we do not, however, have a clear explanation.
Meteorological and glacier mass balance data ( Figure 2 in Sorg et al. (2015), and Bolch and Marchenko (2009)) show that this period covers mostly years with rather cold summer temperatures in the region, and balanced glacier mass budgets, i.e. not 15 particularly warm air temperatures or large melt water productions. Importantly however, and consistent with the glacier mass balances, snow measurements at two stations close to the rock glacier (Bolshoje Almatinskoje Ozero, Great Almatinka Lake , 2516 m asl., and Mynzhilki, 3017 m asl.) (for locations see Bolch and Marchenko, 2009) show strongly above-average snow heights during the 1960s (MeteorologicalService, 1968). Winter temperatures were rather above average in the 1960s (Bolch and Marchenko, 2009). Using a dendrochronological rock glacier activity index, Sorg et al. (2015) show that the period 1960-20 1970 is characterized by particularly high movement activity at Ordzhonikidze rock glacier. Such activity peak is not found for Karakoram rock glacier and another rock glacier in the region studied in the latter work ( Figure S8 in Sorg et al. (2015)).
Increased permafrost and ground ice temperatures due to early snow cover and thus enhanced surface insulation (Bodin et al., 2009), and enhanced amounts of snow melt water (Cicoira et al., 2019b) from exceptionally large snow heights are thus possible explanations for the peak in rock glacier speed during the 1960s. We find no evidence that an earthquake could have 25 caused enhanced displacements (Rosenwinkel et al., 2015). We did not find any indications of particular distortions in the images that could have caused erroneous displacement values. No measurements exist from other studies that could be used to compare our results.
Our time series results for Ordzhonikidze rock glacier are very similar and consistent among all four zones investigated, and most speed variations exceed their error bars. This was much less the case for the photogrammetric results in Sorg et al. (2015), 30 and we suggest our much improved and extended processing as the main reason for the differences, which are in addition for the most part within the large error bars in the 2015 study.
Gorodetsky rock glacier is the only system investigated in this study that is not showing an increase of speed over the study period. Morenny and Archaly rock glaciers show such an increase only in the recent decade and not much exceeding the measurement error bars (Figs 3b and 4a-b). These three rock glaciers are also the ones for which data availability was least good and for which no images existed between 1971 and 2012. Archaly and Morenny rock glaciers show a similar behaviour over time, and increased speed over 1964-1971, similar to Ordzhonikidze rock glacier, which is just 15 km to the northeast, at the northern slope of the Ile Alatau, and likely under similar meteorological conditions. The agreement between Archaly, Morenny and Ordzhonikidze rock glaciers with respect to the increased speeds in the second half of the 1960s, and the high 5 dendrochronological activity during this time on Ordzhonikidze rock glacier (Sorg et al., 2015) suggest that this feature in speed variations could actually be a reliable measurement, possibly related to early snow cover, large snow heights and associated increased ground temperatures and/or amounts of snow melt water. Gorbunov et al. (1992) investigated displacements of boulders on Archaly and Morenny rock glaciers over 1969-1984 and find an overall increase in speeds between 1969-1979and 1979-1984. Over shorter time scales (1973-1977-1978-1982, the authors 10 first identify a decrease in speeds and then again an increase, leading them to suggest the existence of pulse-like speed variations. Our measurements contain no details over these study periods but it seems very reasonable to assume that speeds did actually vary within our 1971-2012 measurement period. (Sorg et al., 2015) did not study the Gorodetsky, Archaly and Morenny rock glaciers. Measurements for some boulders by Gorbunov et al. (1992) on Gorodetsky rock glacier for  reveal speeds similar to our measurements, specifically slightly higher than our 2012-2016 speeds. 15 In-situ measurements of the lower outline of the Gorodetzky rock glacier front revealed that advance rates for  around the location where also we measured the advance are in the order of 0.8 m/a (Gorbunov and Titkov, 1989;Kokarev et al., 1997), i.e. equal to our coarse measurements of 0.8 m/a for 2012-2018, and slightly higher than our 1953-2012 advance of 0.7 m/a. Advance rates increased to about 1.6 m/a in the 1990s (Bolch and Marchenko, 2009). As our 2012-2018 estimates are significantly lower than the latter value for the 1990s, they might be an indication for a decrease in advance rate. Our 1953-20 2012 estimate of 0.7 m/a is, however, not compatible with this assumption unless there were particularly low advance rates between the 1990s and 2012. Given the coarse character of our image-based advance estimates and the high variability of rock glacier advance over short spatial scales, it is difficult to compare our estimates in detail to the previous measurements, though.

Wider geomorphic implications
All photogrammetric velocity fields obtained reflect coherent viscous streams of creeping ice-rich perennially frozen debris. 25 Given the increase in air temperatures in the region during the study period (Bolch, 2007;Sorg et al., 2012;Sorg et al., 2015), the behaviour observed for the bodies of frozen debris in the study area is aligned with findings from Europe and from theoretical considerations that creep speeds of ice-rich perennially frozen debris increase overall under atmospheric warming (Kääb et al., 2007;Delaloye et al., 2010;PERMOS, 2016;Eriksen et al., 2018). Surface speeds for most rock glaciers photogrammetrically investigated in our study were roughly two to four times higher in recent years than during the 1950s and 30 1960s. This indicates a significant increase in sediment and ice fluxes through these rock glaciers, and potentially for the entire study region if one assumes that speed increases apply in overall also to the other rock glaciers or bodies of frozen debris mapped using radar interferometry. A detailed analysis of sediment fluxes is not the focus of this study, but assuming for instance an average thickness of 20 m for the creeping layer (Cicoira et al., 2020) and block movement (i.e. basal velocity equals surface velocity), gives an average increase in volume flux in the order of +200 m 3 /a 2 per rock glacier for the rock glaciers studied in detail. For around 50 large and fast rock glaciers in the study region the latter number would scale to +10'000 m 3 /a 2 . Note, that this is a very rough order of magnitude estimate with a number of unknowns for individual rock glaciers (e.g., from vertical velocity profile and melting ice content) and for the region (e.g., time series representative for all rock glaciers?). 5 Assuming the shrinkage of glaciers in the region (Narama et al., 2010;Farinotti et al., 2015;Severskiy et al., 2016;Treichler et al., 2019) leads to a reduced glacier motion and rate of sediment transport by the glaciers, the increase in rock glacier speed implies that associated sediment transport may currently gain importance relative to sediment transport by glaciers (c.f., Monnier and Kinnard, 2017;Seligman et al., 2019). Such growing geomorphodynamic importance of rock glaciers relative to glaciers under atmospheric warming is in line with the higher thermal inertia of rock glaciers compared to glaciers through 10 their blocky active layer that in essence reduces ground temperatures locally.
The ratios between rock glacier advance rates and surface speeds over the rock glaciers studied in detail cover almost the full spectrum between very small ratios and 1, with most ratios between 0.3 and 0.7. The ratio between both quantities is a function of the ice content that melts out at the front, and the vertical profile of horizontal speeds, i.e. the total material flux into the front (Kääb and Reichmuth, 2005). Unless the ratios are very small (i.e. pointing to very high ice content and/or a shallow 15 moving layer; Morenny), or close to 1 (i.e. low ice content and block movement; Karakoram), we cannot disentangle the influences from ice content and vertical velocity profile on the observed advance rates from our measurements. Relative to each other, Gorodetzky and Ordzonikidze rock glaciers (ratios of 0.7 and 0.6) seem to have a lower ice content and/or less steep decrease of speed down from surface to depth at their frontal parts than Archaly and Kugalan Tash rock glaciers (0.3).

Conclusions 20
Using satellite radar interferometry, we have compiled a slope movement inventory including rock glaciers and other periglacial landforms for the Ile Alatau and Kungöy Ala-Too ranges of northern Tien Shan. We confirm that these ranges host a large number of rock glaciers and identify a large number of comparably fast-moving (>100 cm/a) individual rock glaciers.
The radar-interferometric interpretations derived in this study agree well with photogrammetric offset tracking measurements based on repeat high-resolution optical satellite data, both in terms of overall motion class and surface velocity patterns. When 25 using 12-day Sentinel-1 interferograms, speeds of up to roughly 100 cm/a can be measured, above which phase decorrelation makes precise quantification of speeds impossible. However, the outlines can be drawn using the decorrelation pattern and velocity can be assumed as minimum velocity of e.g. 100 cm/a for a 12-day Sentinel-1 interferogram.
Overall, we find an increase in surface speed for selected rock glaciers. This is in line with theoretical considerations of overall acceleration of creeping ice-rich perennially frozen debris under atmospheric warming, and observational evidence that was 30 so far mainly available for the European Alps. Compared to the dozens of related studies available for the European Alps and the much better availability of suitable data there, our single study can by necessity not provide similarly robust results as available for the Alps, but still represents one of the longest rock glacier speed time-series at regional scale outside Europe.
Our finding suggests that for many rock glaciers of our movement inventory speeds might have been lower during recent decades than contained in this contemporary map. The increase in surface speeds for most time series constructed implies a substantial increase in sediment and ice fluxes through these rock glaciers. This could indicate that sediment transport through rock glaciers and other bodies of frozen debris might currently gain importance compared to a declining sediment transport by 5 glaciers that shrink in the region.
Besides recent acceleration we also find signs of earlier high speeds or increases in speeds . Morenny and Archaly rock glaciers had equally high speeds at the beginning of our measurement series than at their end, Gorodetsky rock glacier likely even higher ones. The most striking peak in speeds is estimated for the second half of the 1960's for Ordzhonikidze rock glacier system on the northern slope of the study region in the Ile Alatau. A possible explanation for this feature are particularly large 10 snow heights in the region during the period of concern that could have contributed to enhanced surface insulation during winter, and associated reduced permafrost cooling and ground ice warming, and/or enhanced flow of snow melt water into the rock glaciers.
All rock glaciers studied photogrammetrically exhibit coherent streams of viscous flow, indicating creep of spatially continuous bodies of ice-rich perennially frozen debris. The good agreement between photogrammetric and radar-15 interferometric results is a strong indication that this finding holds also for the other rock glaciers contained in the kinematic inventory. The rock glacier fronts that reach valley floors show a strongly compressive flow regime, and changes or increases, respectively, in speeds further upstream affect them in a complex way, not only involving changes in horizontal speeds but potentially also in the thickness of the creeping layer of the perennial frozen debris. The fronts might respond rather more passive to upstream speed changes than actively produce own changes. 20 Our study presents one of the most comprehensive data sets of rock glacier kinematics and their changes over recent decades outside of the European Alps. These data may serve a range of follow-up works, such as analysis of rock glacier dynamics, or investigation of associated sediment fluxes in relation to other processes of alpine sediment transport. In terms of additional measurements, the combination of elevation change time series over the rock glaciers with the kinematics of this study might in a next study better illuminate the mass balance of these land forms and their response to climatic changes. 25

Code availability
The image matching code used for this study (Correlation Image AnalysiS, CIAS) is available from http://www.mn.uio.no/icemass (Kääb, 2014). The software used for radar interferometry, the GAMMA remote sensing software, is commercial. 30

Data availability
The rock glacier and slope instability inventory compiled for this study from differential SAR interferograms and optical imagery, 1998-2016, Northern Tien Shan, is available from the Permafrost Information System (PerSys, https://apgc.awi.de/group/persys) (PerSys, 2019). Sentinel-1 data used in our study are freely available from the ESA/EC Copernicus Sentinels Scientific Data Hub at https://scihub.copernicus.eu (Copernicus, 2019). ERS-1/2 and ALOS-1 PALSAR-5 1 data are available from ESA through registration procedures. ALOS-2 PALSAR-2 data are available from JAXA, and TanDEM-X data from DLR through proposal application procedures. Planet data are not openly available as Planet is a commercial company. However, scientific access schemes to these data exist (https://www.planet.com/markets/education-andresearch/). Data from DigitalGlobe satellites (GeoEye, Ikonos, WorldView, Quickbird) and Pléiades are commercial, but programs to facilitate academic access exist. 10

Author contribution
A.Kä. performed the photogrammetric analyses and wrote the paper. T.S. and R.C. performed the radar interferometric analyses and wrote the paper. T.B. assisted with analyses and wrote the paper. H.T. assisted with photogrammetric analyses.
M.S. and A.Ko. contributed to discussions and edited the paper. 15

Competing interests
All authors declare that they have no competing interests.

Acknowledgements
Thanks are due to A. Cicoira, P. Schöneich and L. Rouyet for their careful review comments and W. Haeberli for additional 20 helpful comments that improved the manuscript. We are grateful to the providers of free data for this study: European Space Agency (ESA) / European Commission (EC) Copernicus for Sentinel-1 data, ESA for ERS-1/2 and ALOS-1 PALSAR-1 data, the German Aerospace Center (DLR) for TanDEM-X data, the Japan Aerospace Exploration Agency (JAXA) for ALOS-2 PALSAR-2 data, and Planet for their cubesat data via Planet's Ambassadors Program. ERS-1/2 and ALOS-1 PALSAR-1 data were provided by ESA, courtesy of C1F.6504. ALOS-2 PALSAR-2 images were provided by JAXA, courtesy of RA4-1058. 25 TanDEM-X data was provided by DLR, courtesy of Wegmulle_NTI_INSA3397. We are grateful to CNES/Airbus DS for the provision of the Pléiades satellite data within the Pléiades Glacier Observatory facilitated by Etienne Berthier (LEGOS) and Delphine Fontannaz (CNES).