Climate warming has led to the degradation of permafrost stability in the past half century over the Qinghai-Tibet Plateau

Temperature increases cause a unique type of damage to permafrost. This damage is often expressed as the 10 degradation of permafrost thermal stability, which is very important for engineering design, resource development, and environmental protection in cold regions. This study evaluates the degradation of permafrost stability over the QTP from the 1960s to the 2000s using estimated decadal mean annual air temperatures (MAATs) by integrating remote sensing-based estimates of mean annual land surface temperatures (MASTs), leaf area index (LAI) and fractional snow cover values, and decadal mean MAATs taken at 152 weather stations using geographically weighted regression (GWR). The results reflect a 15 continuous rise of approximately 0.04°C/a in the decadal mean MAAT values over the past half century. Climate warming has led to a reduction in permafrost stability in the past half century. The total degraded area of stability is approximately 153.7610 4 km 2 , which corresponds to 87.98% of the permafrost area in the 1960s. The stability of 75.24% of the extremely stable permafrost, 89.56% of the stable permafrost, 90.3% of the sub-stable permafrost, 92.31% of the transitional permafrost, and 32.8% of the unstable permafrost has been reduced to lower levels of stability. Approximately 49.4% of the 20 unstable permafrost and 95.95% of the extremely unstable permafrost has degraded to seasonally frozen ground. The sensitivity of the permafrost to climate is dependent on its stability level. The mean elevations of the extremely stable, stable, sub-stable, transitional, unstable, and extremely unstable permafrost areas increased by 88 m, 97 m, 155 m, 185 m, 161 m and 250 m, respectively. The degradation mainly occurred from the 1960s to the 1970s and from the 1990s to the 2000s. This degradation has led to increases in risks to infrastructure, increased flood risks, reductions in ecosystem resilience, and 25 positive climate feedback effects. It therefore affects the well-being of millions of people and sustainable development at the Third Pole.

the uncertainty of the gridded air temperatures is significant because of the heterogeneity of the surface characteristics, including snow cover and vegetation, and the locations of weather stations (Vancutsem et al., 2010). Fortunately, the remote sensing era has led to changes in this situation. Thermal infrared remote sensing provides direct observations of land surface temperatures (LSTs) at high spatial and temporal resolutions. For example, the Moderate Resolution Imaging Spectroradiometer (MODIS) LST product is freely available and has been validated over large areas via a series of field 5 campaigns. Its accuracy is better than 1 ℃ (0.5 ℃ in most cases) (Wan et al, 2002;2004;2008). Remote sensing-based estimates of LSTs provide a key high-resolution temperature pattern of the land surface that can potentially be used in monitoring degradation; however, criteria for using LST estimates to distinguish permafrost types are not traditionally available, and the relatively short time series of LST data does not meet the needs of long-term permafrost monitoring. Of the three commonly used predictors for permafrost, the MAGT and permafrost thickness are the most direct indicators of the 10 thermal stability of permafrost (Cheng, 1984); however, long-term measurements of the MAGT and permafrost thickness are almost impossible due to the high cost of drilling boreholes. The MAAT is frequently used in mapping the distribution of permafrost. It is easy to measure and has high spatial representativeness. Importantly, long-term in situ measurements of MAATs are available, and it is possible to estimate MAATs over the QTP using remote sensing-based methods. Therefore, the remote sensing-based LST values can be converted to MAATs and used to monitor the permafrost thermal state, 15 although the potential problems of the MAAT model in predicting permafrost degradation are well known. For example, the performance of the MAAT model is generally affected by the thermal inertia of deep soil layers and geothermal flux (Smith and Riseborough, 2002;Jin et al., 2006;Wu et al., 2010a).
Several previous studies have demonstrated the potential of satellite-based methods in estimating near-surface air temperatures (Hachem et al., 2009;Vancutsem et al., 2010;Yao and Zhang, 2013). The variation in the uncertainty is 20 mainly related to the underlying surface type such as snow cover and vegetation, the amount of solar radiation, and cloud cover (Henderson-Sellers and Hughes, 1982;Zhang, 2005;Vancutsem et al., 2010;Lawrence et al., 2011;Hachem et al., 2012;Ran et al., 2015). Additionally, the highly accurate remote sensing-based snow cover and vegetation products are also available. All of these remote sensing-based data products are very important for estimating the MAAT, which is an air temperature index used in monitoring the thermal stability of permafrost. 25 Therefore, the objective of this study is to evaluate the degradation of permafrost stability over the QTP from 1960 to 2010 by integrating multi-criterion remote sensing observations and an air temperature observation network.

Methods and Datasets
In this study, the degradation of permafrost stability is evaluated based on the MAAT model and the improved MAAT data over the QTP in the past half century. The MAAT in situ measurement data at 152 sites over the QTP and remote sensing 30 data with six independent variables were combined using a Geographically Weighted Regression (GWR) model to estimate the MAAT with a 1 km resolution over the QTP during the past five decades.

Permafrost thermal stability classification system
We use the thermal stability permafrost classification system proposed by Cheng (1984). In this system, permafrost is classified into extremely stable, stable, sub-stable, transitional, unstable, and extremely unstable types, as shown in Table 1.
This system is more useful to describe permafrost degradation from an engineering perspective, rather than changes in the extent of permafrost. The MAAT criterion is available in this system to assess the stability types. On the QTP, a MAAT of 5 −2 ℃ has typically been used to distinguish permafrost from seasonally frozen ground (Cheng, 1984;Ran and Li, 2016). The permafrost stability system was proposed based on the MAGT measurement as an index by analysis of the three-dimensional zonation of the high-elevation permafrost (vertical, latitudinal, and aridity). It is a high level summary of high altitude permafrost zonation. The MAAT index was given according to the statistical relation between MAGT, elevation, and the in situ MAAT measurement (Cheng, 1984). However, the extremely unstable type in the thermal stability 10 classification system proposed by Cheng (1984) refers to regions that include cave ice and frozen gravel below the lower limit of permafrost, which is a very scattered distribution. In this paper, a MAAT of −1 ℃ is simply used to distinguish extremely unstable permafrost from seasonally frozen ground.

Simulation of MAAT using geographically weighted regression
In this study, geographically weighted regression (GWR) is used to simulate MAATs. Local parameters are employed in the 15 GWR model to estimate MAATs while considering the spatial locations of meteorological stations (Brunsdon et al., 1998;Kumar et al., 2012). The weighting is a function of the distance between the location of each regression point and the sites where observations are available. The GWR model used in the present study is shown below in Equation (1): where i y is the MAAT at pixel i, ik x is the k th explanatory factor at pixel i, which is a diagonal matrix; Y is a vector ( 1 n  ) for the dependent variables, i.e., the decadal mean MAAT in the 1960s, 1970s, 1980s, 1990s, and 2000s; and n is the number of MAAT observation stations for each year.
In this study, the Gaussian function is used as a spatial weighting function, as shown in Equation (3): where i d is the distance between the ith observation station and the point to be estimated, and r is the bandwidth parameter.
To accommodate different station densities, the corrected Akaike information criterion (AICc) is used to determine the optimal bandwidth parameters.
A stepwise linear regression analysis is used to select the independent variables for the GWR model ( Table 2). The analysis shows that the use of the MAST, the leaf area index (LAI), the fractional snow cover (FSC), elevation, latitude, and 10 longitude as independent variables results in the highest degree of explanatory power for the past five decades, and the significance level is less than 0.0001. The variance inflation factor (VIF) is used to assess the multicollinearity of the model.
A VIF value <1.5 shows that the degree of tolerance is high, and the multicollinearity of the model is thus acceptable. The performance of the GWR model in the 2010s is shown in Table 2. The five GWR models are then used to estimate the decadal mean MAAT over the QTP for the past five decades. The SAGA (System for Automated Geoscientific Analyses) 15 (Conrad et al., 2015) is used to implement the GWR. Specifically, the GWR for multiple predictor grids geoprocessing tool is used. The Gaussian weighting function and the global search range are used.
Due to the unavailability of the vegetation, snow cover, and LST datasets during the 1960s to 2000s, the effect of the dynamics of vegetation, snow cover, and LST on MAAT during this period is unknown. This will inevitably cause some errors in the estimation of MAAT. Recent studies show that vegetation is increasing overall during the past 30 years, and the 20 snow cover is decreasing overall during the past 15 years over the QTP Huang et al., 2017). The effect of vegetation and snow cover change on MAAT and its feedback process is highly complex. For example, the vegetation-snow interaction effect on MAAT is related to humidity (Zhong et al., 2010;Wu et al., 2015;Yuan et al., 2017).
However, we believe this error mainly occurs at the local level in the nature vegetation dominated areas where the change has occurred at the local level within the last 50 years Huang et al., 2017), and it can be partially 25 compensated by the in situ MAAT measurement over the QTP for the past 50 years.

Evaluation of the degradation of permafrost thermal stability
The following linear regression model is used to evaluate the warming rate or degradation rate over the QTP in the past 50 years: where Y denotes the MAAT or permafrost area, x is the time,  is the error, a is the intercept, and b is the slope (i.e., the warming rate or degradation rate). The statistical significance of the warming rate or degradation rate is evaluated using  Thirteen elevation ranges (<3600 m, 3600-3800 m, 3800-4000 m, 4000-4200 m, 4200-4400 m, 4400-4600 m, 4600-4800 m, 4800-5000 m, 5000-5200 m, 5200-5400 m, 5400-5600 m, 5600-5800 m, and >5800 m) are used to evaluate the elevation dependence of the warming rate. The degradation of permafrost thermal stability is evaluated from two perspectives, the change in area of the different thermal stability types and the spatial heterogeneity of the change. For the area change, we calculate the total area of each thermal stability type during the past five decades and the rate of change (i.e., the degradation 10 rate) using the linear regression model shown in Equation (4). The spatial pattern of the degradation of permafrost thermal stability is evaluated at two levels. At the pixel level, the spatial distribution of the degradation is evaluated. At the level of the thermal stability types, a transfer matrix is used to evaluate the conversions among the thermal stability types (Stehman, 1997). We also analyse the changes in the elevation histograms for each thermal stability type in the past 50 years.

Mean annual land surface temperature
MODIS Terra/Aqua daytime and nighttime LST products (MOD11A1 and MYD11A1, version 5) with a spatial resolution of 1 km and covering 2006 to 2010 were acquired from the Distributed Active Archive Center (DAAC) operated by the U.S. National Aeronautics and Space Administration (NASA). These data are used in this study to estimate MASTs. A pragmatic approach proposed by Ran et al. (2015; is employed to estimate the MASTs using 20 the four daily MODIS LST products. This approach assumes that the arithmetic average of the daytime and nighttime LSTs represent the daily mean LST with acceptable accuracy, and the daily amplitude of LST is more homogeneous than the LST itself (Liu et al., 2006;Kogan et al., 2011;Ran et al., 2015). The approach allows the full use of every value at any time in any pixel of the MODIS LST products through the use of the temporally and spatially complete LST daily amplitude, which is interpolated using a gap filling algorithm (Garcia, 2010). This algorithm employs a penalized least 25 squares regression based on discrete cosine transforms that explicitly utilize information from a time series to predict the missing values. The penalized least squares regression is a thin-plate spline smoother for a generally one-dimensional data array, and it can trade off fidelity to the data versus the roughness of the mean function (Garcia, 2010;Wang, et al., 2012). This approach is easy to implement and independent of other observations. Validation shows that the scheme is effective in restoring the missing values in MODIS instantaneous LST observations and produces a spatially and 30 temporally continuous daily average LST dataset that displays good agreement with observations made at the ground 7 surface. The errors in the results originate mainly from the original instantaneous LST MODIS products. A more detailed description of this scheme can be found in Ran et al. (2015;. The temporally and spatially continuous daily mean LSTs from January 1, 2006 to December 31, 2010 and the corresponding MASTs used in this study are produced using the above approach.

Fractional snow cover 5
Arithmetic mean values of daily cloud-removed FSC products from 2006 to 2010 are used in this study. This product is derived from the daily MODIS 500-m snow cover product (MOD10A1) using a gap filling process based on a cubic spline interpolation algorithm. A comparison with reference FSC data obtained from Landsat Enhanced Thematic Mapper Plus (ETM+) shows the high accuracy with which this product reflects snow cover information over the QTP (Tang et al., 2013).
The cloud-removed FSC products were acquired from the Cold and Arid Regions Science Data Center in Lanzhou, China 10 (http://westdc.westgis.ac.cn).

Leaf area index
Annual mean LAI values obtained from the Global Land Surface Satellite (GLASS), which make up a high-quality LAI product with an eight-day temporal resolution and a 1-km spatial resolution and cover the period from 2006 to 2010, are used in this study. The GLASS LAI product is derived from the fused MODIS and CYCLOPES LAI products, and the 15 remaining effects of cloud contamination have been removed using MODIS time series surface reflectance data and general regression neural networks . The results of validation show that the GLASS LAI product has a lower uncertainty than the MODIS and CYCLOPES LAI products (Xiang et al., 2014). The GLASS LAI product was acquired from the GLASS project website (http://glass-product.bnu.edu.cn/en).

In situ MAAT observations 20
The MAAT measurements, which were collected at 131 stations for the 1960s and 1970s, 133 stations for the 1980s, 144 stations for the 1990s, and 152 stations for the 2000s within the QTP and the surrounding area, were acquired from the data centre of the China Meteorological Administration (http://cdc.nmic.cn). The distribution of the 152 stations for the 2000s is shown in Figure 1. The density of stations in the eastern QTP is higher than other years. The decadal mean MAAT values over the past five decades are used in this study. 25

Validation data
Validation of the long-term stability of permafrost is difficult due to the limited amounts of reference data that are available.
In this study, we evaluate the results by comparing the estimated permafrost distribution in the 2000s with previous regional-scale permafrost maps and borehole measurements at individual sites. The permafrost maps that cover the QTP from , Nan et al. (2002), and Zou et al., (2016) are used at the regional scale. In particular, the map of Zou et al., (2016) integrates the MODIS eight-day LST product within the framework of the temperature at the top of the permafrost (TTOP) model (Smith and Riseborough, 1996), and careful validation of this map has been performed using MAGT data. At the site scale, the MAGT values used in this study were collected from 142 boreholes presented in the existing literature (Yu et al., 2008;Luo et al., 2013)

Ancillary data
The distribution of water bodies in the MODIS land cover product (MOD12Q1) and the map showing the distribution of 10 glacier ice from the second Chinese glacier inventory are used to support the permafrost area statistics. The MOD12Q1 product is used for consistency with the other remote sensing products employed in this study. On the other hand, the glacier extents from the second Chinese glacier inventory are compiled based on Landsat TM/ETM+ images acquired from 2004 to 2011, as well as other ancillary data, such as digital elevation models (DEMs). The robust band ratio segmentation method is first used to delineate the glacier outlines, and intensive manual improvements are then performed to improve its accuracy. 15 An error assessment shows that the area error for all of the glaciers in China is approximately 3.2% (Guo et al., 2015).

Results
Decadal mean MAAT estimates with a resolution of 1 km over the QTP in the past 50 years are produced using the GWR model. The mean coefficient of determination of this model is approximately 0.95. The permafrost stability map in the past five decades is then produced based on the simulated MAAT and the permafrost stability types defined in Table 1. 20

Change of MAAT over the QTP in the past 50 years
The MAATs over the QTP have risen continuously in the past 50 years. The mean MAAT values for the 1960s, 1970s, 1980s, 1990s, and 2000s are -2.38 ℃, -1.85 ℃, -1.78 ℃, -1.32 ℃, and -0.58 ℃, respectively. These values reflect a continuous rise with a rate of approximately 0.04 °C /a. This value is higher than the global average warming rate, as well as the estimated warming rates for the QTP reported by Cheng et al. (2012) and Ran et al. (2016) that are based on interpolated 25 elevation-based air temperature data or surface air temperature reanalysis data. The warming rate in the western part of the QTP is higher than that in the eastern part and depends on elevation, as shown in Figures 2 and 3. The warming rate increases with increasing elevation from approximately 0.33 ℃ per decade at 3600 m to 0.49 ℃ per decade at 5200 m. This finding is similar to that of previous studies (Liu and Chen, 2000;Qin et al., 2009). The physical mechanisms of this phenomenon may be related to the combined effects of the cloud-radiation and snow-albedo feedbacks (Giorgi et al., 1997;Liu et al., 2009;Pepin et al., 2015). These elevated warming rates may have a substantial impact on the thermal stability of the permafrost.

Thermal stability degradation
Based on the map of permafrost stability types covering the past five decades (Figure 4a-e), we analyse the degradation from three perspectives, including temporal changes, spatial changes in the map plane, and spatial changes with elevation. 5

Temporal dynamics of thermal stability
The permafrost thermal stability has degraded continuously over the past 50 years. The area occupied by the stable types has decreased continuously, and the area occupied by the unstable types has increased continuously ( Table 3). The areas occupied by the extremely stable, stable, sub-stable, and transitional types display net decreases of approximately 8.9910 4 km 2 (72.79%), 27.0610 4 km 2 (70.12%), 9.3010 4 km 2 (27.24%), and 1.1810 4 km 2 (4.77%) from the 1960s to the 2000s, 10 respectively. In particular, the stable type displays the most serious degradation, and its rate of loss is approximately 6.1510 4 km 2 (15.94%) per decade. Moreover, the area occupied by the unstable type has increased by approximately 3.9910 4 km 2 (9.02%) at a rate of 1.0610 4 km 2 (2.4%) per decade. Specifically, this degradation mainly occurred during the 1960s to 1970s and the 1980s to 1990s for the extremely stable type, the 1960s to the 1970s and the 1990s to the 2000s for the stable type, and the 1980s to the 2000s for the sub-stable type. The area occupied by the extremely unstable type has not changed 15 substantially. Overall, the warming climate has caused a degradation of permafrost stability. If glaciers and the extremely unstable type are included, the total permafrost area has decreased significantly from 174.7610 4 km 2 in the 1960s to 133.110 4 km 2 in the 2000s at a rate of approximately 9.5210 4 km 2 (5.45%) per decade, and this loss of area occurred mainly during the 1960s to the 1970s and the 1990s to the 2000s (Table 3).

Spatial changes in thermal stability 20
The degradation of thermal stability has occurred over a broad region of permafrost on the QTP within the past 50 years, especially during the 1960s to the 1970s and the 1990s to the 2000s. The degradation of permafrost stability in the western QTP was serious during the 1960s to the 1970s. In the subsequent 40 years, the degradation of permafrost stability in the QTP was relatively homogeneous (Figure 4f-i). Specifically, the extents of the extremely stable, stable, and sub-stable types retreated from the south to the north (Figure 4a-e). The extents of the transitional, unstable, and extremely unstable types 25 extended northward correspondingly. Approximately 42.30% of the area occupied by the extremely stable type, 42.09% of the area occupied by the stable type, and 39.83% of the area occupied by the sub-stable type have degraded to the stable, sub-stable, and transitional types from the 1960s to the 1970s, respectively. At the same time, approximately 57.26% of the area occupied by the transitional type, 29.34% of the area occupied by the unstable type, and 59.47% of the area occupied by the extremely unstable type, have degraded to the unstable type, the extremely unstable type, and seasonally frozen ground, 10 respectively. Overall, approximately 75.24% of the area occupied by the extremely stable type, 89.56% of the area occupied by the stable type, 90.3% of the area occupied by the sub-stable type, 92.31% of the area occupied by the transitional type, and 32.8% of the area occupied by the unstable type have degraded to lower levels of stability in the past 50 years (Table 4).
The reduction in the area occupied by the permafrost is mainly due to the degradation of the area occupied by the unstable and extremely unstable types. Approximately 49.4% of the area occupied by the unstable type and 95.95% of the area occupied by 5 the extremely unstable type has degraded to seasonally frozen ground (Table 4). The total degraded area is approximately 153.7610 4 km 2 , which accounts for 87.98% of the area occupied by the permafrost region in the 1960s (Figure 4j). The area of permafrost for which the stability has not changed is approximately 2110 4 km 2 (12.02%). This area is mainly distributed in the central part of the plateau, which contains extremely high mountains, and it is dominated by the extremely stable type.
Notably, the stability of a permafrost area of approximately 1.6310 4 km 2 has increased. This area is found primarily east of 10 Lhasa in the southeastern part of the QTP, which is a major centre of marine glaciers and snow cover in China (Figure 4j).
The increased permafrost stability in this area may have large uncertainties; the uncertain MAAT trend is estimated using regression parameters that are appropriate for low-elevation areas, due to the lack of long-term MAAT measurements in the high mountain regions where glaciers and snow are prevalent. The effects of snow or glacier cover may be more important than those of the MAAT. Although records of long-term snow cover and glacier changes in the past 50 years are not 15 available in this study, the sensitivity of glacier and snow cover in a warming climate is dependent on the climate zone. Low snow-climate sensitivities have been found in continental interior climates with relatively cold and dry winters (Brown and Mote, 2009). Larger glaciers have lower climate sensitivities (Ding and Haeberli, 1996;Ye et al., 2001).
Additionally, the complex process and limited knowledge for permafrost-glacier interactions may enhance the uncertainty (Haeberli, 2005;Otto and Keuschnig, 2014). Therefore, we believe the permafrost stability in this area has not changed 20 substantially in the past 50 years, based on this low climate sensitivity.

Elevation changes in permafrost stability type distributions
The elevation statistics of the distribution of the permafrost stability types over the QTP in the past five decades indicate that the elevation occupied by each permafrost stability type in the QTP has increased continuously (Table 5). For the extremely stable type, the mean elevation of the distribution decreased from 5240 m to 5161 m from the 1960s to the 1970s and then 25 rose continuously at a rate of approximately 56.4 m per decade. The reduction in elevation is mainly due to the degradation of the extremely stable permafrost type in the Kailas Mountains. This caused the fluctuation of the mean elevation for extremely stable permafrost during the 1970s to 1980s and reduced its statistical significance (low R in Table 5) for the increasing rate of mean elevation over the past 50 years. Overall, in the past 50 years, the mean rate of increase of the extremely stable type has been approximately 24.7 m per decade. Moreover, the mean elevation of the stable, sub-stable, 30 transitional, unstable, and extremely unstable types have risen at a rate of 23.6 m, 36.3 m, 43 m, 36.5 m, and 56.2 m, respectively. Overall, the mean elevation of the extremely stable, stable, sub-stable, transitional, unstable, and extremely unstable types increased by 88 m, 97 m, 155 m, 185 m, 161 m, and 250 m, respectively, over the past 50 years. This result indicates that the climate sensitivity of permafrost is dependent on the stability level. The extremely unstable permafrost type is the most sensitive of the permafrost types to climate warming. As in the last section, the degradation mainly occurred from the 1960s to the 1970s and from the 1990s to the 2000s.

Cross validation and uncertainty analysis 5
We validate the permafrost extent only in the 2000s because long-term records of permafrost stability and extent are not available in earlier periods, as mentioned in section 2.4.5. Comparison of the estimated permafrost extent in the 2000s with the permafrost map provided by Zou et al., (2016) shows that the difference is small. Within permafrost areas, the extremely unstable type of permafrost mainly refers to cave ice and frozen gravel, which are distributed below the lower limit of permafrost (Cheng, 1984). This kind of permafrost is usually not counted in the total area of permafrost. Therefore, the 10 permafrost area in the 2000s is approximately 107.19 10 4 km 2 if glaciers and lakes are neglected. This result is similar to that of Zou et al., (2016), who showed that the permafrost area in the 2000s was approximately 106.4710 4 km 2 . The permafrost distribution is also very similar to that presented by Zou et al., (2016) (Figure 5b). The consistency between the two distributions is 92%, and the kappa coefficient is approximately 0.82. At the site scale, 89% of the 142 locations are consistent with the borehole survey, whereas this proportion is only 74%, 28%, and 86% for the maps of Cheng, (1996), 15 Nan et al. (2002), and Zou et al. (2016), respectively. These proportions indicate that the accuracy of the permafrost extent identified in this study is at least comparable with that of Zou et al. (2016).
The uncertainty of the results may result primarily from the MAAT model, insufficient resolution, inaccuracies in the surface station data, or the sparseness of these stations, which are especially sparse in high mountain areas. First, the response time and the depth to which permafrost is affected by climate warming depend on the extents, durations, amplitudes, and rates of 20 climate warming and are closely related to soil types, surface coverage, ice content, groundwater occurrence, geothermal anomalies, and human activities (Stieglitz et al., 2003;Zhang, 2005;Lawrence et al., 2008;Cheng and Jin, 2013;Westermann et al., 2016). For example, the low heat conductivity of soil leads to lags between increases in surface temperatures and the subsequent increases in permafrost temperature or reductions in permafrost thickness . The delay time is longer for permafrost thickness than temperature and varies with the thermal stability type Wu et al., 2010a). 25 For the stable type, the degradation of permafrost may be delayed by "thermal offset" and "seasonal offset" effects in the permafrost table due to the negative heat budget; i.e., the amount of heat released from the active layer during the winter is greater than the amount of heat absorbed in summer (Smith and Riseborough, 2002;Wu et al., 2010a). For the unstable type, a positive heat budget appears in the upper soil layer that leads to a greater degradation rate than that seen in stable permafrost, since the thickness of the unstable type is smaller than that of the stable type Wu et al., 2010a). However, the 30 complex physical mechanisms of the interactions between climate change and permafrost are currently poorly understood (Jin et al., 2011), and a large degree of uncertainty may exist in previous evaluations as well as the permafrost area change over the past 50 years in this study. Despite current warming, large permafrost areas may persist due to the thermal inertia of permafrost (Cheng et al., 2012). Second, the thawing of the base of the permafrost induced by the geothermal heat flux leads to the permafrost degrading from bottom to top (Jin et al., 2006;Wu et al., 2010a). The MAAT model cannot reflect the change of geothermal flux from the crustal interior. Additionally, the geothermal flux data are generally limited or unavailable. The missing geothermal heat flux may lead to a delay in permafrost degradation, especially for the stable 5 permafrost, because the geothermal flux is independent of air temperature. Third, although the resolution of the simulation has been significantly improved to 1 km, it is still coarse relative to the degradation rate of mountain permafrost. The degradation of mountain permafrost is presented in terms of the increase in the elevation of the lower limit of the permafrost, which is generally approximately a hundred metres. A 1 km change in the horizontal extent change may correspond to a change in elevation of hundreds of metres. Last, the lack of long-term MAAT measurements in the glacier-and 10 snow-dominated high mountain regions may lead to errors in the estimated MAATs.
Overall, the "accelerated degradation" effect of the MAAT model may be partly counteracted by the "delayed degradation" effect of the missed geothermal heat flux. Long-term observation shows that the mean increasing rate of ground temperature at 10-20 m depth in the QTP is approximately 0.024 ℃ (Zhao et al., 2010;Wu et al., 2010b;Jin et al., 2011), which is comparable with the warming rate of air temperature. This shows that the evaluation results of permafrost stability 15 degradation using the MAAT model is generally accepted at the overall QTP scale.

The implications of the degradation of thermal stability
The degradation of permafrost stability in the QTP has important impacts on the safety of infrastructure in permafrost regions, water quality, ecosystem health, and the feedbacks on regional and global climates. First, as the permafrost stability degrades, the risk of deterioration and damage to engineered structures in permafrost zones will increase. This indicates that the 20 measures used to prevent permafrost degradation may need to be enhanced for new structures. For example, permafrost accounted for 90.1% of a 10-km-long segment of the QTR from Golmud to Lhasa in the 1960s, and these permafrost areas were dominated by the sub-stable type; however, after 50 years (i.e., in the 2000s), these permafrost areas accounted for only 67.77% and were dominated by the unstable type. For "warm" permafrost areas that are dominated by unstable permafrost, an enhanced measure to prevent permafrost degradation, i.e., the proactive roadbed cooling approach, has been successfully 25 applied in constructing the QTR (Cheng, 2004;Cheng et al., 2008). Second, the degradation of permafrost in the QTP may affect the hydrologic cycle in the Third Pole region, which includes the QTP and the surrounding arid regions.
Permafrost controls the distribution, recharge, flow paths, discharge, dynamics, and hydrochemistry of groundwater (Cheng and Jin, 2013). The degradation of permafrost affects the interactions among the surface water, subsoil water, and groundwater by changing the hydraulic conductivity and hydraulic connectivity of the soil. The degradation of the ice-rich 30 permafrost itself makes important contributions to surface runoff and the development of thermokarst lakes in the inner Tibetan Plateau . The enhanced drainage may lead to increases in flood risk (Larsen et al., 2008) and reductions in ecosystem resilience via seasonal shifts in stream flow and groundwater abundance, because the decrease in permafrost water storage capacity in the QTP will lead to a reduction in dry-season water availability. All of these changes will affect the well-being of millions of people and sustainable development at the Third Pole, which contains the headwater areas of several of the major rivers in Southeast Asia, such as the Yellow, Yangtze, Mekong, Yarlung Zangbo and Shiquan Rivers.
The Third Pole also includes many inland rivers, such as the Shiyang, Heihe, Shule, and Tarim Rivers, in northwestern China.
Last, the permafrost region in the QTP contains approximately 160 Pg of organic carbon (Mu et al., 2015) and many 5 thermokarst lakes and wetlands (Niu et al., 2011;Luo et al., 2015). Thawing of the permafrost may lead to the disappearance or growth of thermokarst lakes (Smith et al., 2005), which may further affect greenhouse gas emissions and produce a feedback effect on climate change (Tarnocai et al., 2009;Schuur et al. 2009;Schaefer et al., 2011;McCalley et al., 2014).

Conclusions
This study evaluates the stability degradation of permafrost over the QTP from the 1960s to the 2000s based on the estimated decadal means of the mean annual air temperatures (MAATs) over the Qinghai-Tibet Plateau (QTP) in the past 50 years obtained by integrating remote sensing-based mean annual land surface temperatures (MASTs), leaf area index (LAI) and fractional snow cover values, and decadal mean MAATs measured at 152 weather stations using a geographically weighted 15 regression (GWR) model. Cross validation shows that the accuracy of the estimated permafrost extent is greater than that of previous maps.
The decadal mean MAATs reflect a continuous rise at a rate of approximately 0.04 °C /a during the past half century. The warming rate increases with increasing elevation from approximately 0.33 ℃ per decade at 3600 m to 0.49 ℃ per decade at 5200 m and then decreases as elevation increases further. Climate warming has led to the degradation of permafrost stability in 20 the past half century. The area occupied by the stable permafrost types has continuously decreased, and the area occupied by the unstable permafrost types has continuously increased. The total degraded area is approximately 153.7610 4 km 2 , which accounts for 87.98% of the permafrost area in the 1960s. The stability for all permafrost types have degraded to lower levels.
The extent of the extremely stable, stable, and sub-stable types retreated from the south to the north, whereas the extent of the transitional, unstable, and extremely unstable types extended northward. The mean elevations of the extremely stable, stable, 25 sub-stable, transitional, unstable, and extremely unstable types increased by 88 m, 97 m, 155 m, 185 m, 161 m and 250 m, respectively. This result indicates that the climate sensitivity of permafrost is dependent on the stability level. The degradation mainly occurred during two periods that include the 1960s to the 1970s and the 1990s to the 2000s. The degradation of permafrost stability in the QTP has important impacts on the safety of infrastructure, flood risks, ecosystem resilience, and climate feedbacks, as well as the well-being of millions of people and sustainable development at the Third The uncertainties inherent in this analysis cannot be discounted. These uncertainties are due to asynchronous changes in near-surface air temperatures and deep soil layer temperatures, the missing geothermal flux, insufficient resolution, or the inaccuracies and sparseness of the surface station data employed. As this evaluation is empirically based, obtaining more convincing results requires additional data, especially from the deep layers of soils. The development of new, fast, and inexpensive sensors and robust machine learning methods will assist in this effort. A physically based definition of 5 permafrost stability and an improved physically based model will contribute to the prediction of permafrost stability degradation and its interactions with the engineering stability of infrastructure, the water cycle, and climate change.
Table 1. Classification system used to assess permafrost stability (Cheng, 1984) Type Mean annual ground temperature (℃)