Monitoring Regional-Scale Surface Deformation of the Continuous Permafrost in the Qinghai–Tibet Plateau with Time-Series InSAR Analysis

Abstract

As an important indicator of permafrost degradation, surface deformation is often used to monitor the thawing and freezing process in the permafrost active layer. However, due to the large area of the continuous permafrost of the Qinghai–Tibet Plateau (QTP) and the large amount of data processed by conventional time-series InSAR, previous studies have mostly focused on local area investigations, and regional characteristics of surface deformation of the continuous permafrost area on the QTP are still unclear. In this paper, we characterized surface deformation in space and time over the main continuous permafrost area on the QTP, by analyzing 11 ascending and 8 descending orbits of Sentinel-1 SAR data acquired between 2018 and 2021 with the time-series InSAR processing system LiCSAR. The reliability of the InSAR deformation results was verified by a combination of leveling measurement data, the intercomparison of overlapping area results, and field verification. The results show that the permafrost regions of the central QTP exhibited the most significant linear subsidence trend. The subsidence trend of permafrost on the QTP was mainly related to the thermal stability of permafrost, and the regions with larger subsidence rates were concentrated in sub-stable, transitional and unstable permafrost areas. We also found that, according to analysis of time-series displacement, the beginning and ending times of permafrost thawing were highly spatially heterogeneous, with the time of maximum thawing depth varying between mid-October and mid-November, which was probably attributed to the active layer thickness (ALT), water content in the active layer, and vegetation cover in these regions. This study is of great significance for understanding the changing trend of permafrost on the QTP under the background of climate change. In addition, this study also demonstrates that combination of Sentinel-1 SAR images with the LiCSAR system has significant potential for detecting permafrost deformation with high accuracy and high efficiency at regional and global scales.

1. Introduction

The Qinghai–Tibet Plateau (QTP), the third pole of the Earth, contains the largest permafrost area to be found in middle and low latitudes of the world [1]. As the water tower of Asia, the QTP has a large amount of ground ice stored in its permafrost area, with a reserve of about 12.7 trillion cubic meters [2]. In addition, the permafrost region of the QTP is also a huge organic carbon reservoir, and the organic carbon storage within 2 m of the surface of the permafrost region is about 17 Gt [3]. Under the background of global warming, the QTP is an amplifier of global warming. In the past 40 years, the rate of temperature increase of the QTP has been twice that of the global rate during the same period, resulting in the severe degradation of permafrost [4,5]. The degradation of permafrost on the QTP has an important impact on regional and global water and carbon cycles [6,7].

The permafrost degradation is mainly manifested in a reduction in the permafrost distribution area, an increase in ground temperature, the thinning of permafrost thickness, the thickening of the active layer, and the melting of underground ice [8,9]. The active layer is the buffer layer between the permafrost and the ground surface, which thaws in summer and freezes in winter, with seasonal changes, making the ground surface subside in summer and uplift in winter [10,11]. In addition, the melting of shallow ground ice in permafrost areas causes a linear subsidence trend on the surface, which is particularly obvious in the ice-rich permafrost area [12]. Therefore, monitoring the surface deformation of permafrost areas is vital to investigate permafrost degradation [13].

Time-series InSAR has proven an established remote sensing technique to detect small surface deformations with high precision, and at high resolution, and has been successfully applied to monitor surface deformation in the permafrost regions of Siberia, Alaska, and the QTP [14,15,16,17,18,19,20]. In the Alaska North Slope region, Liu et al. used ERS data to obtain the surface deformation results based on time-series InSAR technology, indicating that time-series InSAR technology can provide a new means of observation for permafrost system dynamics and permafrost state monitoring [14]. Chen et al. used time-series InSAR to observe the freeze–thaw cycle process of the active layer of the Alaska North Slope, and found that the deformation observed by InSAR was mainly related to the soil moisture content and the freeze–thaw cycle of the active layer [15]. In the Siberia region, Chen et al. used Sentinel-1 data to obtain the results of surface deformation in the Yedoma region, indicating that there was a significant seasonal thawing subsidence in summer, and they observed that there is greater seasonal thawing subsidence in flat areas than in sloped areas [16]. Abe et al. used L-band InSAR to monitor surface subsidence due to thermokarst in central Yakutia, and the results showed that InSAR has good application for monitoring thermokarst processes in permafrost areas [17].

In the permafrost region of the QTP, Daout et al. used the time-series InSAR method to obtain surface deformation results for the northwestern QTP during 2003–2011 and separated its long-term deformation rate and seasonal deformation amplitude, finding that the areas with a large seasonal deformation amplitude are mainly concentrated in the sedimentary basin area [18]. Li et al. used the time-series InSAR technology to obtain time-series deformation results for the Wudaoliang area of the Qinghai–Tibet Railway and verified the InSAR deformation results by using the GPS leveling data of 24 points. The results show that the time-series InSAR technology can be effectively applied to the surface deformation monitoring of permafrost in the QTP [19]. Chen et al. used Sentinel-1 data to obtain the time-series deformation results of the area along the Qinghai–Tibet Railway, and combined them with environmental factors such as terrain, vegetation, and ground ice to analyze the main influencing factors of linear deformation rate and seasonal deformation amplitude in the study area [20]. On the QTP, current research on the deformation of permafrost is mainly carried out in the Qiangtang Basin or along the Qinghai–Tibet Railway. However, there is a lack of research on the deformation of large-scale permafrost areas in the complex environment of the QTP.

Recently, the number of satellites with InSAR capability and the volume of associated data have increased significantly, and the improvements in data quality and in processing methods have greatly improved the capability and application range of InSAR in ground deformation monitoring [18,19,20]. In particular, the Sentinel-1 satellite constellation is capable of conducting earth observations with a revisit period of 6–12 days, making it an important data source in ground deformation monitoring [21,22]. In addition, the recently developed open-source InSAR data-processing platform LiCSAR system, with the help of extensive Sentinel-1 data, can more conveniently obtain InSAR deformation results, which provides a valuable platform for large-scale InSAR deformation research [22,23].

This study aimed to characterize regional-scale surface deformation in space and time over the main continuous permafrost area of the QTP by using Sentinel-1 SAR datasets and the LiCSAR time-series InSAR processing system. A total of 11 ascending and 8 descending orbits of Sentinel-1 SAR data acquired between 2018 and 2021 were exploited to detect the large-scale surface deformation of the continuous permafrost area of the QTP. Furthermore, based on the surface deformation results validated by levelling measurements of 21 surveying points, the spatial distribution characteristics of permafrost deformation were analyzed. Finally, the potential of these time-series permafrost deformation measurements for possible application to large-scale permafrost degradation investigations on the QTP was considered.

2. Study Area

In this paper, the main continuous permafrost area of the QTP is taken as the study area. The permafrost region of the QTP is the highest-elevation and largest permafrost region in the middle and low latitudes of the world [24]. The elevation of the study area is mainly concentrated between 4000 and 5500 m, the terrain is relatively flat, and the slope is restricted to between 0 and 25°. The climate of the study area is characterized by long sunshine hours, strong solar radiation, low precipitation and low temperature throughout the year [25]. The vegetation types in the study area are mainly alpine meadows, alpine grasslands and alpine deserts [25].

The continuous permafrost on the QTP is mainly distributed around the Qiangtang Basin, and there is also a continuous permafrost area in the southern part of the Qilian Mountains [1]. The permafrost area on the QTP is about 1.06 × 10⁶ square kilometers, accounting for 40.2% of the total area (Figure 1b) [26]. According to the mean annual ground temperature (MAGT) in the permafrost region, the extremely stable permafrost (MAGT < −5 °C) is mainly distributed in the Karakoram and Altun Mountains; the stable type (−5 °C < MAGT < −3 °C) and the sub-stable type (−3 °C < MAGT < −1.5 °C) are mainly distributed in Hoh Xil, Kunlun Mountain Pass and other areas; the transitional type (−1.5 °C < MAGT < −0.5 °C) and unstable type (−0.5 °C < MAGT < 0.5 °C) are mainly distributed in the east and south of the permafrost area of the QTP (Figure 1a); and the proportions of transitional and unstable permafrost are 36.77% and 20.69%, respectively [27]. The active layer thickness on the QTP is mainly between 0.8 and 3.5 m, and gradually increases around the Qiangtang Basin. The ALT in the Qiangtang Basin is between 0.87 and 1.78 m, with an average value of 1.02 m, while in the Beiluhe Basin, the ALT is about 1.5–3.5 m [28,29]. Affected by the continuous increase in temperature on the QTP, the active layer has been thickening and accelerating in recent years, and the permafrost has been degraded significantly [1].

3. Datasets and Methods

3.1. SAR Data

The LiCSAR system is an open-source InSAR data processing system released by the Centre for Observation and Modelling of Earthquakes, Volcanoes and Tectonics (COMET), which is mainly used for large-scale interferometry processing based on Sentinel-1 data [23]. The LiCSAR system can automatically provide geocoded and unwrapped interferograms and coherence maps with a resolution of 0.001 degrees (WGS-84 coordinate system). The differential interferometry data provided by the LiCSAR system platform are processed by InSAR data-processing software GAMMA and Snaphu. GAMMA software is used for data preprocessing and differential interferometry processing, and Snaphu is used for unwrapping processing. Before the unwrapping processing, pixels with coherence less than 0.5 are rejected [23].

Data for each frame are obtained after differential interferometry based on Sentinel-1 data, including unwrapping differential interferometry images, coherence images and DEM images of the corresponding area. In order to obtain the deformation results for the permafrost on the QTP, a total of 19 Sentinel-1 data were used, including 11 ascending-orbit data and 8 descending-orbit data. The data acquisition period was mainly concentrated in 2018–2021. The information of the data used is shown in

Table 1. The Sentinel-1 data information used in this paper.

Frame	Date	Ascending/Descending	Sentinel-1 Data Count
099A_05217_131313	12/01/2018–14/01/2020	Ascending	58
099A_05416_131313	07/01/2019–08/01/2021	Ascending	58
026A_05128_131313	21/01/2017–26/01/2019	Ascending	60
143A_05651_131313	04/03/2018–05/03/2020	Ascending	54
041A_05629_131313	22/12/2018–23/12/2020	Ascending	58
041A_05430_131313	11/10/2018–12/10/2020	Ascending	54
012A_05443_131313	07/03/2018–08/03/2020	Ascending	55
085A_05455_131313	12/03/2018–13/03/2020	Ascending	58
085A_05654_131313	24/03/2018–06/02/2020	Ascending	54
158A_05451_131313	03/05/2018–18/03/2020	Ascending	59
056A_05457_131313	07/10/2019–10/02/2021	Ascending	43
077D_05487_141313	20/02/2017–10/02/2019	Descending	42
077D_05685_121313	07/09/2018–23/05/2020	Descending	44
004D_05513_131313	05/01/2018–07/01/2020	Descending	53
150D_05505_131313	04/04/2019–24/03/2021	Descending	60
121D_05470_131313	09/01/2019–10/01/2021	Descending	53
121D_05668_131313	17/10/2018–02/07/2020	Descending	40
121D_05217_131313	09/01/2019–10/01/2021	Descending	62
019D_05619_131313	02/01/2019–28/12/2020	Descending	48

.

3.2. Leveling Data

In order to verify the validity of the deformation results for the permafrost region obtained by time-series InSAR, we obtained the long-term deformation velocity of 21 points along the Qinghai–Tibet Railway measured by leveling. The location distribution of leveling observation points is shown in Figure 1. We obtained the long-term deformation velocity of embankments at 21 points along the Qinghai–Tibet railway from 2005 to 2011 by using a high-precision leveling instrument. Although the observation period of long-term deformation velocity measured by leveling is inconsistent with that of InSAR, the railway embankments are in a relatively stable state due to various protective measures taken during the construction of the Qinghai–Tibet Railway [30,31]. Therefore, we can still verify the validity of time-series InSAR results by comparing the long-term deformation velocity measured by leveling with the linear deformation velocity calculated by InSAR.

3.3. Time-Series InSAR Processing

In this study, we used Sentinel-1 data and the LiCSAR system for time-series InSAR processing and, on this basis, used the permafrost deformation model to calculate the linear deformation velocity. The data processing flow is shown in Figure 2. The LiCSAR system provides a software package LiCSBAS which can be used for time-series InSAR processing. During time series analysis, statistical quality checks are carried out by calculating the average coherence and the percentage of the number of effective unwrapping pixels in the interferogram to identify and eliminate bad data [32]. The acquired unwrapped images may still contain a large number of unwrapping errors, and if these data are not rejected, they can cause significant errors in the deformation results. According to Equation (1), if there is no unwrapping error in the three interferograms, the loop phase should be close to zero [32]. If all loops associated with an interferogram are problematic, the interferogram concerned is considered to be one that may contain many unwrapping errors, and the interferogram is rejected in the next processing step.

$${\mathsf{\Phi }}_{123}={\mathsf{\Phi }}_{12}+{\mathsf{\Phi }}_{23}-{\mathsf{\Phi }}_{13}$$

(1)

The checked unwrapped correct interferograms were used to calculate the time-series deformations. Although Sentinel-1 data have good temporal and spatial baselines, due to severe de-coherence, periods of time with no acquisitions, and the results of network refinement, gaps are prone to exist in the small baseline network. In order to avoid this problem, this study used the NSBAS method with time constraints on deformation time series for time series analysis [22,32]. After the time series analysis, we obtained the initial time-series deformation results. However, there are also relatively large pixel-to-pixel differences in a stack of unwrapped data, and the deformation results for some pixels are not reliable due to unwrapping errors. To solve this problem, we used the noise indices generated in the previous steps to mask out the “bad” pixels and obtained the time-series deformation results after masking.

In large-scale InSAR deformation processing, space-related tropospheric delay, ionospheric delay, orbit error and other noise signals are the main error sources. We used a spatiotemporal filter (i.e., high pass in time and low pass in space) that can be applied to attempt to separate these noise components from the displacement results [32]. In addition, the deformation signals caused by tectonic movement usually also appear as spatially correlated signals. Following previous study [11], the spatiotemporal filtering can effectively reduce the influence of tectonic deformation signals. In Figure 3a, there is an obvious deformation trend caused by tectonic motion, while Figure 3b shows the deformation result after spatiotemporal filtering, and the tectonic deformation trend has been effectively removed. Finally, the results of time-series deformation caused by permafrost deformation are obtained.

Due to the frost heave and thaw subsidence of the active layer, the deformation of the permafrost presents a periodic trend, and the surface deformation also exhibits a linear subsidence trend with the melting of ground ice near the permafrost table. Therefore, the deformation model of the permafrost region can be expressed as Equation (2), and the linear deformation rate can be calculated by Equation (2) [18,33].

$$\mathrm{LOS}\left(t_i\right)=\mathrm{V}{t}_i+a_1\times \sin \left(\frac{2\pi }{T}{t}_i\right)+a_2\times \cos \left(\frac{2\pi }{T}{t}_i\right)+\mathrm{c}$$

(2)

where $\mathrm{LOS}\left(t_i\right)$ is the deformation in the LOS direction, V is the linear deformation velocity, T is the deformation period (assumed to be one year), $t_i$ is the time span from the starting image, $a_1$ and $a_2$ are the seasonal deformation parameters, and $\mathrm{c}$ is the residual term.

4. Results

4.1. InSAR Permafrost Deformation

According to the results for linear deformation velocity in the permafrost area of the QTP (Figure 4), large linear subsidence mainly occurs in the transitional and unstable permafrost area in the center of the QTP. The linear subsidence trend in this area is generally greater than 10 mm/yr, and even more than 20 mm/yr in some areas. The average ground temperature of permafrost in this area is relatively high, and with the increase in temperature, the degradation of permafrost has become more intense, and the thickness of the active layer shows a significant thickening trend [1]. The statistical model results show that the average change rate of ALT along the Qinghai–Tibet highway from 1981 to 2018 reached 19.5 cm/10 yr [34]. In addition, the ground ice content in this area is large [29], and the more severe degradation of permafrost can lead to the severe melting of ground ice near the permafrost table. Therefore, the large linear subsidence trend in the transitional and unstable permafrost areas in the middle of the QTP is probably caused by the melting of ground ice and sedimentary compaction after the thickening of the active layer.

In the west of the study area, the linear subsidence trend in most areas is around 5 mm/yr, and in some areas between 5 mm/yr and 10 mm/yr. Here the overall deformation trend is relatively stable, and the areas with large linear subsidence trend are distributed exclusively around lakes and rivers. When combined with the distribution results of permafrost stability types in the QTP, we find that the region is mainly composed of sub-stable and transitional permafrost, and the annual average ground temperature is mostly less than −0.5 °C. Compared with the permafrost area in the central QTP, the degree of degradation of permafrost in this area is relatively light, while the thickening rate of the active layer and the melting degree of ground ice are relatively small. Therefore, the overall linear subsidence trend in this area is relatively stable, and the large subsidence trend is obvious only in the ice-rich permafrost area covered by alluvial sediments around lakes and rivers.

In the northeast of the study area, the region is mainly composed of sub-stable-type permafrost (−3 °C < MAGT < −1.5 °C), and the overall linear subsidence trend is relatively small. The larger linear settlement rate mainly occurs in the sedimentary area around Hala Lake and its subsidence rate exceeds 10 mm/yr, which may be caused by the degradation of permafrost around the lake and the compaction of alluvial sediments [11].

We specified a study area with a linear subsidence velocity greater than 10 mm/yr and counted the areas and area ratios of different stable types of permafrost in the study area. The results are shown in

Table 2. Area statistics table of different stability types of permafrost in the area with linear subsidence rate greater than 10 mm/yr.

Permafrost Type	Area (km2)	Proportion
Extremely stable	122	0.3%
Stable	3745	8.9%
Sub-stable	13,732	32.9%
Transitional	16,992	40.7%
Unstable	7202	17.2%

. In the study area, the area with a linear subsidence rate greater than 10 mm/yr is about 41793 km², of which the area of transitional permafrost is about 16992 km², accounting for 40.7% of the total, followed by sub-stable and unstable permafrost regions, which account for 32.9% and 17.2%, respectively. This result shows that the obvious linear subsidence trend in the study area is mainly found in the transitional, sub-stable and unstable permafrost, and this may be related to the severe degradation of permafrost in these areas. This result also shows that the linear subsidence trend of permafrost is related to the thermal stability of permafrost, indicating the importance of permafrost deformation monitoring for the study of permafrost changes under the background of climate change.

Compared with previous local deformation monitoring research [11,18,20], the overall deformation trend of the main continuous permafrost regions of the QTP can now be revealed. We found that the linear subsidence trend in the central region of the QTP is most significant in the main continuous permafrost regions of the QTP. The region is mainly an ice-rich permafrost area, and important transportation routes such as the Qinghai–Tibet Railway and Qinghai–Tibet Highway pass through the region. Therefore, under the background of climate change, this region is one of the most noteworthy and continuously monitored areas in the permafrost region of the QTP.

Due to the large scope of the study area, in order to analyze the deformation of different areas, we selected four points with large linear deformation rates along the east–west direction: P1, P2, P3 and P4 (Figure 4). Figure 5 shows the results of local linear subsidence velocities at P1, P2, P3 and P4. Among them, the three regions of P1, P2, and P3 are all transitional permafrost, and the P4 region is dominated by sub-stable-type permafrost. P1 is located in a basin area between Tanggula Mountain and its glaciers, in which thermokarst lakes are densely distributed, and the Qinghai–Tibet Highway passes through the basin (Figure 5a). The linear subsidence velocity in the basin exceeds 20 mm/yr, and the distribution of thermokarst lakes is very dense (Figure 6a). The formation of thermokarst lakes is mainly caused by ground ice melting, which also confirms that the large linear subsidence velocity in this area may be due to the thawing of ice-rich permafrost. The linear subsidence velocity in both the P2 and P3 regions exceeds 5 mm/yr (Figure 5b,c), but the linear subsidence velocity in the P2 region is obviously greater than that in the P3 region, which may be related to the fact that the P2 region is closer to the southern boundary of permafrost and its permafrost degradation is more severe. The linear subsidence in the P4 area mainly occurs in the alluvial sediment area around Hala Lake, most obviously in the the alluvial sediment area to the south, and the linear subsidence rate generally exceeds 10 mm/yr (Figure 5d). There are extensive thermokarst lakes and polygonal cracks in this area, which also reflects the degradation of ice-rich permafrost in this area (Figure 6b). In addition, since the area is composed mainly of alluvial sediments, the larger linear subsidence trend may be related to sediment compaction in the area in addition to the impact of ground ice melting.

4.2. Result Verification

4.2.1. Intercomparison of Overlapping Area Results

In order to quantitatively verify the validity and spatial consistency of the deformation velocity results, a total of six overlapping areas were identified for verification, including four overlapping areas of the same track and two overlapping areas of different tracks (Figure 4). We selected regions with slopes of less than 10° within the overlapping regions for intercomparison to reduce the influence of geometric factors. In all intercomparison areas, the deformation rates obtained from different images are processed by difference [35]. Figure 7 shows the statistical histograms of the difference between different overlapping areas. In the six overlapping areas identified, the mean value of the linear deformation velocity difference is within 1 mm/yr, and the standard deviation is within 5 mm/yr, indicating that the deformation results obtained from different images have good consistency.

4.2.2. InSAR Results vs. Leveling Measurements

The slope of all our leveling verification points is less than 5°, and in Figure 7e, we present the comparison results of ascending and descending tracks in the area of leveling verification point Y13. The deformation results of ascending and descending tracks have good consistency; therefore, we think that the deformation at the leveling verification points is mainly caused by vertical deformation. In order to compare the linear deformation velocity calculated by InSAR with the long-term deformation velocity of leveling measurement, we divided the LOS-direction deformation result by the cosine of the incident angle and converted it into the vertical-direction deformation result [36,37]. The long-term deformation velocity of 21 leveling observation points along the Qinghai–Tibet Railway is compared with the linear deformation velocity calculated by InSAR. The results are shown in

Table 3. Comparison results of the long-term deformation rate measured by leveling and the linear deformation rate of InSAR.

Observation Point	Leveling Rate (mm/yr)	InSAR Rate (mm/yr)	Difference between Leveling and InSAR (mm/yr)	Observation Point	Leveling Rate (mm/yr)	InSAR Rate (mm/yr)	Difference between Leveling and InSAR (mm/yr)
Y1	1	−4	−5	Y12	0	0	0
Y2	−1	−4	−3	Y13	−1	0	1
Y3	0	3	3	Y14	−4	−7	−3
Y4	−7	−11	−4	Y15	−7	−12	−5
Y5	−4	−2	2	Y16	−11	−13	−2
Y6	−5	−8	−3	Y17	−2	−3	−1
Y7	−2	−3	−1	Y18	−2	−9	−7
Y8	4	3	−1	Y19	−4	−6	−2
Y9	−16	−13	3	Y20	−4	−9	−5
Y10	−14	−7	7	Y21	−1	−5	−4
Y11	−6	−8	−2	Mean ± std			3 ± 2

. For all 21 leveling observation points, the absolute value of the difference between the long-term deformation velocity and InSAR linear deformation velocity at 19 monitoring points is less than or equal to 5 mm/yr, with only two monitoring points showing a value of 7 mm/yr. The mean value of the absolute value of the deformation velocity differences for all 21 monitoring points is 3 mm/yr, and the standard deviation is 2 mm/yr. The comparison results show that the deformation results of the time-series InSAR measurement agree well with the leveling results.

Although our leveling verification points are linearly distributed, the overlapping areas we use for intercomparison are uniformly distributed overall. Combining the leveling verification results and the intercomparison results of the overlapping areas, we can still indicate that our deformation results are valid.

4.2.3. Field Verification

In order to deeply analyze and verify the deformation results, we carried out field verification on the QTP. Figure 8 shows the field photos and time-series deformation results for the field investigation regions on the QTP.

Figure 8a shows the field photo of Wudaoliang (WDL) region, and Figure 8b shows the time-series deformation results for this area. Obvious thermokarst collapse is observed in this area, which is mainly formed by the degradation of ice-rich permafrost. Its time-series deformation is obviously periodic and accompanied by a large linear subsidence trend, and the linear subsidence rate reaches 14.9 mm/yr. Figure 8c,e show the field photos of Fenghuoshan (FHS) and Sewuxiang (SWX) regions, respectively. There are retrogressive thaw slumps caused by the degradation of ice-rich permafrost in both areas. Figure 8d,f show the time-series deformation results for FHS and SWX, respectively, and the linear subsidence rates are 9.2 mm/yr and 12.5 mm/yr, respectively. The obvious thermokarst landforms in these field investigation regions confirmed that the large linear subsidence in this area may be caused by the degradation of ice-rich permafrost, and the field investigation results verify the validity of the deformation results.

5. Discussion

5.1. Comparisons with Previous Studies

In the permafrost region of the QTP, many studies have used the time-series InSAR method to monitor surface deformation. To further evaluate the validity of our results, we compared the areas where our results coincided with previous studies, and the results are shown in

Table 4. Table of comparison between the results of this study and previous studies.

Authors	Study Area	SAR Data	Observation Period	Linear Rate (mm/yr)	This Study Result (mm/yr)
Chen et al.	central QTP	Sentinel-1	2014–2019	−15 to 15	−20 to 15
Wang et al.	Beiluhe Basin	TerraSAR	2014–2016	−15 to 0	−15 to 5
Li et al.	Wudaoliang region	ENVISAT	2006–2009	−15 to 15	−20 to 10

. Chen et al. used Sentinel-1 data to obtain deformation results for the central Qinghai–Tibet Plateau from 2014 to 2019 [20]. The overall deformation rate in the region is −15 mm/yr to 15 mm/yr. Our results are similar to Chen’s results for the overall deformation trend of the central QTP, and our deformation rate is −20 mm/yr to 15 mm/yr, which is consistent with their results. Wang et al. used TerraSAR data to obtain deformation results for the Beiluhe Basin from 2014 to 2016 [29] and found that the overall deformation rate in the basin is −15 mm/yr to 0 mm/yr, which is consistent with our results. Li et al. used ENVISAT data to obtain the deformation results for Wudaoliang region from 2006 to 2009 [19] and found that the total deformation rate of this area is between −15 mm/yr and 15 mm/yr, which is consistent with the results of this paper. In general, the results of this study are consistent with the previous studies, and the small differences may be related to the differences in the observation periods of the data, which further illustrates the validity of the results of this paper.

5.2. Spatial Heterogeneity of Time of Maximum Thawing Depth

The active layer exhibits frost heave and thaw subsidence with changes in temperature, i.e.. frost heave in winter and thaw subsidence in summer. On such a large scale as the permafrost area of the QTP, due to the great differences in temperature, topography and ALT in different areas, the characteristics of surface deformation in these areas also differ obviously. We obtained time-series deformation information for four regions—P1, P2, P3 and P4—and the temperature data for the corresponding regions during the deformation observation time (the temperature data are ERA-5 reanalysis data), as shown in Figure 9, and analyzed the differences in surface deformation characteristics in the different regions.

Figure 9a shows the time-series deformation and temperature change curves for the point P1. The daily mean temperature in this area reaches its annual maximum from mid-June to mid-September. The mean temperature during this period is above 6 °C. The mean daily temperature drops sharply after mid-September and reaches minimum in mid-January, after which the temperature begins to rise again. Correspondingly, the deformation time series of point P1 shows the beginning of thaw at the end of April, maximum depth of subsidence reached at the end of October, the beginnings of freezing at the end of November, with maximum uplift reached at the end of February which remained basically until mid-April.

The time-series deformations of the four points were analyzed., Their thawing start times were similar, and they all began to thaw at the end of April (Figure 9). They differed in the times at which they reached their maximum thawing depth. P1 reached the maximum thawing depth between the end of October and mid-November, P2 from mid-October to the end of October, P3 at the end of October, and P4 in mid-October. The time difference between the four points in reaching the maximum thawing depth is 15–30 days.

According to the summer thawing process of the active layer, the thawing process of the active layer is a one-way process, the heat transfer direction of the active layer moves from top to bottom, and the thawing front gradually migrates downward and gradually reaches the maximum thawing depth [38,39,40]. Therefore, the time taken to reach the maximum thawing depth is greatly affected by the ALT. In theory, the deeper the ALT, the longer the thawing time. Therefore, the differences in time taken to reach the maximum thawing depth at the four points may be related to the ALT. However, the ALT is not the only factor affecting time taken to reach the maximum thawing depth. The water content in the active layer and vegetation cover are also important factors affecting time differences in reaching the maximum thawing depth [41,42].

5.3. Application Potential of Large-Scale InSAR Deformation

Against the background of climate change, the degradation of permafrost on the QTP has intensified, leading to changes in the surface deformation characteristics of permafrost. Time-series InSAR can monitor the dynamic process of permafrost freeze–thaw cycles. Therefore, the large-scale time-series InSAR results for permafrost on the QTP have great application potential for the monitoring of permafrost degradation and the inversion of related parameters.

The ALT is an important attribute of permafrost, and the monitoring of ALT can reflect the degradation state of the permafrost. However, the current remote-sensing inversion of the ALT of the QTP mainly uses surface temperature data for inversion, which has a low spatial resolution and has difficulties in reflecting local detailed information. The time-series InSAR deformation results can reveal the seasonal deformation amplitude of the active layer, and the combination of the seasonal deformation amplitude, soil porosity and soil moisture content can be used for the inversion of the active layer thickness [10]. Therefore, the large-scale InSAR deformation results for permafrost in the QTP can reveal the seasonal deformation amplitude of permafrost, and then be used for the high-resolution inversion of ALT on the QTP.

Permafrost degradation leads to the intensification of ground ice melting in ice-rich permafrost, and the melting of ground ice leads to obvious long-term subsidence of the surface. Time-series InSAR results can accurately reflect the linear deformation of permafrost except seasonal deformation and can be used to monitor land subsidence caused by ground ice melting. Therefore, the large-scale InSAR deformation results for the permafrost area of the QTP have important application potential for monitoring the changes of ground ice reserves in the permafrost areas of the QTP.

6. Conclusions

We used Sentinel-1 data and the LiCSAR system to obtain time-series deformation results for the main continuous permafrost areas of the QTP. Comparing the InSAR results with the deformation results of 21 leveling measurement points, the mean of the absolute value of the difference between the InSAR deformation rate and the leveling deformation rate is 3 mm/yr. The results of comprehensive leveling data analysis, intercomparison of overlapping areas and the field investigation show that the InSAR deformation results are reliable. We found that the permafrost regions of the central QTP exhibited the most significant linear subsidence trend. The linear subsidence trend in the permafrost region was mainly related to the thermal stability of permafrost, and more than 90% of the regions with linear subsidence rates greater than 10 mm/yr in the study area were concentrated in the sub-stable, transitional and unstable regions. In addition, we found that according to analysis of time-series displacement, the beginning and ending time of permafrost thawing was highly spatially heterogeneous, reaching maximum thawing depth between mid-October tand mid-November, which was probably attributed to the active layer thickness (ALT), to water content in the active layer, and to vegetation cover in these regions.

Under the background of climate change, the permafrost on the QTP is undergoing a severe degradation process, and this study is an important reference for understanding the changes in the major continuous permafrost regions on the QTP. In addition, the large-scale deformation monitoring results in permafrost regions have important application potential for the inversion of permafrost parameters.